Article Cite This: J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
pubs.acs.org/JCTC
Formulation of Small Test Sets Using Large Test Sets for Efficient Assessment of Quantum Chemistry Methods Bun Chan* Graduate School of Engineering, Nagasaki University, Bunkyo 1-14, Nagasaki 852-8521, Japan
Downloaded via DURHAM UNIV on July 13, 2018 at 20:54:29 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.
S Supporting Information *
ABSTRACT: In the present study, we have examined in detail literature data of deviations for a wide range of (mainly) DFT methods for the extensive MGCDB82 set (∼4400 data points) of main-group thermochemical quantities. We use the data and standard statistical techniques (lasso regularization and forward selection) to devise the MG8 model for linearly combining assessment results of a collection of small data sets to accurately estimate the MAD of MGCDB82. The MG8 model contains a total of 64 data points representing noncovalent interactions, isomerization energies, thermochemical properties, and barrier heights. It is thus well suited for rapid evaluation of new quantum chemistry procedures. We propose that a value of ∼4 kJ mol−1 for an estimated MAD by the MG8 model (EMADMG8) to be an initial indicator of a highly robust quantum chemistry method, with large deviations occurring mainly for properties (such as heats of formation) that are known to be difficult to accurately compute. For methods with larger EMADs, we emphasize the importance of more thorough testing, as these methods are likely to have a larger number of outliers, and it may be less trivial to anticipate circumstances under which large deviations occur. In relation to this aspect, we have applied the same generally applicable statistical techniques to further formulate small-data-set models for assessing the accuracy for some properties that are not covered by MG8 nor by MGCDB82. They include the MOR13 model for metal−organic reactions, the SBG5 model for semiconductor band gaps, and MB13 for stress-testing methods with artificial species.
■
INTRODUCTION The development of low-cost yet adequately accurate quantum chemistry methodologies, notably density functional theory (DFT) methods,1 has come a long way in recent decades. Lately, their formulation relies quite heavily on fitting to and/ or testing with large sets of benchmark data.2−6 During the early days of this trend, moderately sized data sets such as G27,8 have been used, for example, for the optimization of the B97 functional,9 on which many modern DFT methods are based. This set of thermochemical quantities contains roughly 150 data points. The size of typical test sets has grown since then. Very recently, Mardirossian and Head-Gordon have compiled the MGCDB84 data set for main-group chemistry,10 and it astonishingly comprises nearly 5000 data points! Likewise, the recently introduced GMTKN55 set11 contains approximately 1500 thermochemical quantities for main-group molecular species. This is still very substantial even though it is relatively modest when compared with MGCDB84. The availability of these extensive data sets has on the one hand provided developers of DFT and other approximate methods an excellent means for thorough testing and tuning of their methods. On the other hand, when one is at an early stage of formulating a new procedure, it would be advantageous to explore a wide range of methodological directions in a more efficient manner. Thus, a common practice is to use a smaller training set for development purpose and then test the resulting method with a larger test set. However, even the size © XXXX American Chemical Society
of training sets has grown substantially. For instance, the newly developed ωB97M-V6 functional was initially trained with a subset of MGCDB84 that already contains close to 900 data points, before eventually being tested on the full set of nearly 5000 benchmark values. In our own development of new quantum chemistry procedures, including both DFT methods12−16 as well as high-level composite protocols,17−23 we typically use smaller training sets such as the E218,19 and E324−26 sets that we compiled. It is thus of interest to us, and we believe also to other developers, to evaluate whether one can construct models of small data sets to represent larger ones. Indeed, some of the widely used small data sets (e.g., AE627 and DBH2428) are derived from larger parent sets. Whether one can substantially downsize very large sets such as MGCDB84 and GMTKN55 without compromising the utility, however, remains an open question. In the present study, we attempt to use the results related to MGCDB84 to devise smaller and generally applicable data sets, with the hope of facilitating the development of new methodologies in a timely fashion.
■
COMPUTATIONAL DETAILS Standard computational chemistry calculations were carried out with Gaussian 1629 and Q-Chem 5.0.30 Corrections for Received: May 25, 2018
A
DOI: 10.1021/acs.jctc.8b00514 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation dispersion effects were obtained using the dftd3 program.31 Geometries and benchmark relative energies were obtained from previous studies. Notably, the ones in the MGCDB84 set were from ref 10. In general, we used the 6-311+G(3df,2p) basis set32 in our single-point energy calculations. For the calculation of species in data sets that contain transition metals, the def2-TZVP basis set33 was used for all systems in that data set. For periodic systems, we used the m6-311G* basis set.34 The MGCDB84 data set contains 84 subsets with a total of nearly 5000 data points. Among these, two subsets (AE18 and RG10) are generally excluded when assessing accuracies of the various DFT methods because they are more suitable for the training of methods but their relevance to practical applications is rather limited.10 In the present study, we adopt this convention and do not include them in our analysis. To avoid confusion, we will refer to this data set (with approximately 4400 data points) as the MGCDB82 set. Statistical analyses were carried out with JMP Pro 14.35,36 Unless otherwise noted, energies in the text are given in kJ mol−1.
actions on the basis of number of data points. However, the large magnitudes of the thermochemical relative energies (REs), in conjunction with the still fairly large numbers of data points, make them the major contributor on the energy basis. They are also associated with large variations in relative energies, as indicated by the large standard deviations (SDs), both in absolute terms and relative to their RE values. In the extensive benchmark study with MGCDB82, mean absolute deviations (MADs) for the eight categories of maingroup thermochemical data have been determined, but no overall gauge was presented.10 We can obtain an overall MAD by summing the corresponding weighted MADs, with the coefficient for each group determined by the number of data points that it contains. This is equivalent to averaging the absolute deviations for all data points. While it is quite common to use such a straightforward metric for assessing a method, alternative measures are employed in other cases (root-mean-square deviation, geometric mean, etc.; see for example refs 2 and 37). We have considered and compared some of the other metrics and come to the conclusion that the choice of using one or the other may not be of critical importance (see the Supporting Information). In the present study, we use MAD as we prefer direct inspection of the actual magnitudes of the deviations. Guided Derivation of Smaller Data Set Using Statistical Model-Building Tools. Let us now begin our discussion on building a model with a substantially smaller data set than MGCDB82. Our general tactic, irrespective to the exact strategy employed, is to treat the overall MAD for MGCDB82 as a linear combination of the MADs for the 82 subsets (plus a constant term for better quantitative agreement). An exact relationship exists when all subsets are included, but it becomes approximate when some are excluded. Our initial target is to eliminate as many data points as possible by surveying the importance of the subsets as contributors to MGCDB82. In the process of reducing the number of subsets, we aim for the model of smaller sets to achieve two objectives. First, it should reproduce as closely as possible the MADs for the collection of DFTs covered in ref 10. Second, we aim to retain some information about the characteristic strengths and weaknesses of the DFT methods, though not at the same level of detail as what one would obtain from MGCDB82 itself. In the first instance, we use, as much as possible, standard statistical model-building tools to select the important subsets. This approach has an appealing feature in that it reduces the requirement for judicial decisions and is therefore seemingly more objective than some of the strategies that we will discuss later. In our preliminary study, we have examined lasso and elastic net regularizations, as well as forward selection and backward elimination regressions using a variety of program options.38 Our discussion will focus on the results obtained with lasso and forward selection techniques, as the general features of their selection outcomes are representative of those obtained with regularization and stepwise methodologies. Attempts in Subset Selection Based on Lasso Regularization and Forward Selection. Our first attempt is to use the lasso selection technique to drastically reduce the number of subsets while maintaining a good correlation between the estimated MADs (EMADs) and the actual MADs for MGCDB82 (MADMGCDB82) using the coefficient of determination (R2) as the indicator. Remarkably, it requires just a combination of three subsets to reach an R2 that is larger than 0.98 (Table 2, Model 1). The use of nine subsets will take
■
RESULTS AND DISCUSSION Description of the Existing MGCDB82 Data Set That Will Form the Basis of New Models of Smaller Test Sets. Before we devise models of new and smaller data sets for assessing computational chemistry methods, we will take a step back and provide a brief description of the data that we utilize for our analysis. Specifically, in the present study, we will focus mainly on the MGCDB82 set (see computational details above for its definition) due to its inclusion of a large amount of results and raw data. A perspective for the range of systems contained in MGCDB82 is provided in Table 1. The 82 Table 1. Composition of the MGCDB82 Data Set Including Brief Descriptions of Its Eight Categories, the Numbers of Data Points (#), Average Magnitudes of the Reference Relative Energies (RE), and Standard Deviations Between the Data (SD)a data set MGCDB82 NCED NCEC NCD IE ID TCE TCD BH
description noncovalent interaction (easy, dimer) noncovalent interaction (easy, cluster) noncovalent interaction (difficult) isomerization energy (easy) isomerization energy (difficult) thermochemistry (easy) thermochemistry (difficult) barrier height
#
RE
SD
4399 1744
187.4 19.3
724.0 47.4
243
207.3
266.1
91
41.7
54.1
755 155 947 258 206
24.0 106.8 573.7 534.3 100.4
18.5 158.7 1247.3 1664.5 100.0
Energies are given in kJ mol−1.
a
subsets of MGCDB82 have been categorized into eight groups.10 Among these eight, three represent noncovalent interactions, two are compendiums of isomerization energies, another two contain collections of thermochemical quantities, and one represents barrier heights. We can see that almost half of the data points (#) are related to noncovalent interactions. The remaining half is dominated by isomerization energies and thermochemistry, with the latter contributing a larger proportion. The number of barrier height data constitutes less than 5% of the total. One may perceive that MGCDB82 is over-represented by noncovalent interB
DOI: 10.1021/acs.jctc.8b00514 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
of building a model for predicting the overall MAD for MGCDB82 involves successive forward selection regressions. We start by carrying out an unconstrained regression with all 82 subsets, and we trace the resulting solution path from a large collection of subsets toward smaller compendia. Along this path, we locate the point at which the elimination of a subset is associated with the elimination of the category that it belongs to. This subset is protected from being eliminated in our next (constrained) forward selection. We iterate this procedure until we reach a collection of eight subsets each representing one unique category. This model is shown in Table 3, along with some of its characteristics.
Table 2. Some Models Produced by Lasso Regression for Estimating the Overall MADs of MGCDB82 as Linear Combinations of Its Subsets for the ∼200 Computational Quantum Chemistry Methods Assessed in Ref 10, Groups and the Subsets That They Contain, and the Associated R2 Values for the Fits model
group
data setsa
R2
1
NCED TCE NCED NCEC IE TCE TCD BH NCEC
S66x8 HAT707nonMR, TAE140nonMR S66x8 NCEShields38, H2O20Bind10 YMPJ519 AlkAtom19, HAT707nonMR, TAE140nonMR HAT707MR HTBH38 FmH2O10, WATER27, Shields38, 3B-69TRIM, H2O20Bind10 NBC10, S66x8, Ionic43 Bauza30 YMPJ519 Styrene45, C20C24 AlkAtom19, G21EA, G21IP, BDE99nonMR, HAT707nonMR, TAE140nonMR, BSR36, HNBrBDE18 HAT707MR, TAE140MR, PlatonicTAE6, PlatonicIG6, PlatonicHD6 NHTBH38, HTBH38, BHPERI26, CR20, CRBH20, PX13WCPT27
0.984
2
3
NCED NCD IE ID TCE TCD BH
0.998
Table 3. Data Set Devised Using Forward Selection for Estimating the Overall MADs of MGCDB82, the R2 Values for the Correlations with the Overall MADs of Their Parent Categories, and Number of Data Points
1.000
a
group
data seta
R2
#
NCEC NCED NCD IE ID TCE TCD BH
Shields38 AlkBind12 CT20 H2O20Rel4 Styrene45 TAE140nonMR PlatonicIG6 NHTBH38
0.918 0.951 0.143 0.507 0.711 0.985 0.125 0.865
38 12 20 4 45 124 6 38
See ref 10 and references therein for the definition of these data sets.
See ref 10 and references therein for the definition of these data sets.
the R2 value to 0.998 (Model 2). However, this model has a drawback. As mentioned above, the subsets in MGCDB82 are classified into eight categories (Table 1). We deem such a classification quite reasonable and the inclusion of data representing each group a favorable feature of MGCDB82. Unfortunately, the nine-subset model (Model 2 in Table 2) contains only six of the eight groups. We have briefly examined the use of different options for lasso regression. A comparison of several resulting solution paths shows that, to include at least one subset from each category, the smallest collection would contain ∼30 subsets (e.g., Model 3). In our opinion, they still represent a considerable hindrance to effective initial screening/assessment of theoretical methods. The results in Table 2 show that the variations in the accuracies for MGCDB82 for the various theoretical methods can be well described by the variations in the accuracies for just a few subsets. The lasso technique, being a “strong eliminator”, successfully identified this characteristic and provides us with a simple model but at the expense of a lessened diversity in the categories of chemical properties represented by MGCDB82. In this case, we deem the lasso technique (and even the milder eliminator, elastic net) being not optimal for helping us to retain some minor subsets. This is because a lasso solution path tends to involve large-step changes, with many minor predictors being eliminated at the same time. This also eliminates the opportunity to rank and select some of these. Next, we will discuss the use of the alternative forward selection technique to gradually build a small model with all categories in MGCDB82 being included. A characteristic of the forward selection method is that each step along the solution path involves addition (or deletion, in the reverse direction) of just one factor, which facilitates straightforward fine-tuning of the resulting model. Our method
The subsets included in this model contain 287 data points, which is substantially less than the full MGCDB82 set. The overall EMADs for MGCDB82 for all methods examined correlate well with the actual MADs, with an R2 value of 0.989. At this point, let us reiterate that we wish to formulate a model of small data sets that can (1) estimate MADMGCDB82 values with good accuracy and (2) yield information about the characteristics of the assessed methods. The second target forms the basis of including a distinct data set to represent each of the eight categories in MGCDB82. A relevant question regarding this aspect is whether the eight sets are representative of the categories to which they belong. A survey of the relevant R2 values shows that, while some sets in this model correspond quite well with their parent groups, this is not the case for others. In particular, the data sets for the NCD and TCD categories are associated with poor R2 values of 0.143 (NCD) and 0.125 (TCD). The data set that represents the IE group is also not very representative of its parent group (R2 = 0.507). We therefore seek to formulate an alternative model in which the subsets that it contains show better correlations with their respective categories. Formulation of a Predictive Model with Representative Subsets for the Different Categories. To select what we deem the most representative subset for each category, we rank all the subsets in the same group using their R2 values for the correlations between their MADs and those for their parent groups using all the assessed quantum chemistry methods. The eight sets chosen by this means are shown in Table 4. They are S66x8 (NCEC), 39 H2O20Bind10 (NCED), 40 Bauza30 (NCD),41 YMPJ519 (IE),42 C20C24 (ID),43 TAE140nonMR (TCE),44 TAE140MR (TCD),52 and DBH24 (BH).28 Among this collection of eight subsets, five of them have R2 values that are larger than 0.95, and two have R2 values that are between 0.90 and 0.95. For the remaining one, its R2 value is
a
C
DOI: 10.1021/acs.jctc.8b00514 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation Table 4. Data Sets with the Largest R2 Values for the Correlations with the Overall MADs of Their Parent Categories and Number of Points in the Original (#) and Trimmed (#*) Sets group
data seta
R2
#
#*
NCEC NCED NCD IE ID TCE TCD BH
S66x8 H2O20Bind10 Bauza30 YMPJ519 C20C24 TAE140nonMR TAE140MR DBH24
0.997 0.974 0.913 0.973 0.746 0.985 0.958 0.917
528 10 30 519 8 124 16 24
14 3 9 18 4 5 4 7
We believe that we have at this point reached our target. The data set size in our model is fairly small. It estimates the overall MADs of MGCDB82 with good accuracy. The correlation is shown pictorially in Figure 1, with the standard error being 0.9
a
See ref 10 and references therein for the definition of these data sets.
0.745. This last R2 value is by no means excellent, but it is still quite adequate as a qualitative/semiquantitative measure for the accuracy of a theoretical method for computing this category of properties. If we were to represent each category with better accuracy using models of multiple subsets, our analysis shows that, with criteria of R2 > 0.995 and 0.990, we would need totals of 30 and 24 subsets, respectively. One could of course loosen the criterion to further reduce the number of subsets, but the choice is somewhat arbitrary, and we decide to use the minimal set of eight for its simplicity and uniqueness. We linearly combine the MADs of these eight sets to yield overall EMADs for MGCDB82. The coefficients and the constant term are determined by a straightforward linear regression. We can also derive alternative fitting terms by considering the number of data points represented by each category. A comparison of the two sets of parameters is given in the Supporting Information. Here we simply note that the results thus obtained are not very different, and we use the regression values for the somewhat better fit to the target values. The formula for estimating the overall MADs for MGCDB82 is given by eq 1, and the R2 for the correlation between the estimated and actual MADs is 0.992.
Figure 1. Mean absolute deviations (MADs, kJ mol−1) for MGCDB82 versus those estimated with the (64-data-point) MG8 model for a set of 200 quantum chemistry methods assessed in ref 10.
kJ mol−1. In addition, it yields reasonable information about how a computational chemistry method performs for the eight categories of properties contained in the MGCDB82 set. These outcomes are accomplished mainly by statistical optimization, with (in our opinion) few prejudiced decisions being involved. We believe that the approaches and results presented in this work may be useful for the construction and improvement of data sets more generally. We designate our optimal 64-datapoint model as MG8 and the MADs estimated by this model as EMADMG8. We have provided the fitted parameters in a spreadsheet within the set of Supporting Information files, and the same spreadsheet also enables automated calculation of the EMADs for the overall MGCDB82 set and its subsets. Cross-Validating MG8 with Alternative VM6 and VM9 Models. Let us now consider the matter of further probing the applicability of MG8. To determine the quality of EMADMG8 as a proxy for MADMGCDB82, an ultimate test would be to devise and/or assess a series of unrelated methods first with MG8 and then compare the EMADs with the actual ones. This approach necessitates a large number of computations using the entire MGCDB82 set. Instead, we propose an approximate strategy as a means to cross-validate the use of MG8 as a benchmark for new and existing quantum chemistry methods. In this approach, we devise other small models that are independent of MG8 but are also based closely on MGCDB82. Thus, prospective agreements between the results obtained with these models would lend further support to the generality of MG8. In this case, as we intend for these alternative models to simply be sanity-checks for MG8 rather than ones to be used more generally, we statistically choose the subsets without setting any criterion on their types. In the first instance, we select a model with the smallest number of data sets that provides an R2 value above a relatively more modest target of 0.990. Because we aim for the validation model to be independent of MG8, the eight data sets contained in the MG8 model are excluded from our regression. The regression leads to six data sets. For each of these sets, we reduce the number of data points in a similar manner as being done for MG8, i.e., by generalized regressions with a criterion of an R2
EMAD = 0.4533 × MADS66x8 + 0.0078 × MADH2O20Bind10 + 0.0580 × MADBauza30 + 0.2084 × MAD YMPJ519 + 0.0116 × MADC20C24 + 0.0720 × MADTAE140nonMR + 0.0232 × MADTAE140MR + 0.0278 × MADDBH24 + 2.52 (kJ mol−1)
(1)
Although this model represents MGCDB82 and its eight categories reasonably well, a disadvantage is that it contains a significant number of data points (1259) because some of the subsets are quite large. To tackle this challenge, we ask the following question: if a collection of 82 data sets can be represented by just eight subsets, can the overall characteristic of a large set be described by a much smaller number of data points? We thus seek a representation of the MAD for each set as a linear combination of the absolute deviations of its individual data, with a form that is similar to eq 1 above. We examine the use of lasso and forward selection to downsize the eight sets and choose regression solutions with the least number of points while satisfying the condition of having an R2 that is larger than 0.995. In this manner, the number of data points in the entire model is reduced to just 64, and the EMADs still correlate well with the actual ones with an R2 of 0.988. D
DOI: 10.1021/acs.jctc.8b00514 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
Table 5. Estimated MADs for MGCDB82 by MG8 (with an Uncertainty of ±1.8 kJ mol−1), VM6 (± 2.2) and VM9 (± 0.8) for Four Existing (HSEB-D to APF-D) and Two New (MBh-D and PGTh) DFT Methods and the Final Estimated Values Based on These Three Models MG8 VM6 VM9 all
HSEB-D
EDF2-D
WP04-D
APF-D
MBh-D
PGTh
5.9 7.0 6.3 6.6 ± 2.5
9.8 10.1 8.3 9.8 ± 2.4
12.1 10.3 13.3 11.2 ± 3.0
8.3 10.4 8.6 9.5 ± 3.0
6.6 6.9 6.5 6.9 ± 2.1
14.2 12.9 14.7 13.4 ± 2.6
5). They suggest a moderate increase in general accuracy in the order PGTh → WP04-D → EDF2-D ∼ AFP-D → MBh-D ∼ HSEB-D. Importantly, although the three independent models contain in total less than 200 data points, they are drawn from 23 subsets of MGCDB82, which is a small but not insignificant portion of the total collection of 82 subsets. Overall, we deem the cross-validations given by VM6 and VM9 supportive of the MG8 model as a means for estimating the likely MAD for the much larger MGCDB82 set. By extension, it provides as a useful indicator for gauging the general accuracy of a DFT method for main-group molecular chemistry. As an aside, we note that a slightly revised HSE method (HSE-HJS)49 has an MAD of 8.8 kJ mol−1 for MGCDB82,10 and this somewhat larger MAD when compared with that for HSEB-D is consistent with our previous findings.15 Beyond MG8, Lessons Learned from the Full MGCDB82 Set. In previous sections, we have illustrated the utility of using the small test sets in the MG8 model for probing the accuracy of quantum chemistry methods. Significantly, the results are by-design indicative of what one might expect from the much larger MGCDB82 set. With the small model providing an efficient and fairly robust means for assessing the general accuracy of a method, one might be tempted to bypass further assessments in the formulation or selection of a method. While this may not be unreasonable in some cases, one must bear in mind that the general accuracy (overall and for specific categories) of a method is not a completely conclusive indicator for a particular problem. In this section, we will use the data given in ref 10 to briefly evaluate how one might cautiously approach such an issue. Among the (∼200) methods assessed in ref 10, ωB97M-V6 has the lowest overall MAD of 3.3 kJ mol−1, whereas the functional ranked 100th is B97-D50 (MAD = 9.3 kJ mol−1), which we consider a method representing an “average” DFT. Let us now examine the distributions of their deviations for the MGCDB82 set (Figure 2). For ωB97M-V, there are ∼80% of
being larger than 0.990. The resulting model contains just 34 data points. This model is designated as validation model 6 (VM6). The R2 value for the correlation between the estimated and actual MADs for MGCDB82 is 0.982, which is not as good as the one for MG8 but still quite reasonable. The standard error in the estimation is 1.1 kJ mol−1. We have also formulated another model in which the subsets in MG8 and VM6 are being excluded in the regressions, but a stricter criterion of 0.998 for R2 is used. The resulting model contains nine subsets. Upon trimming the size of each subset, the final model consists of 84 data points. This second validation model is termed VM9. It has an R2 of 0.997 for the correlation with MADMGCDB82 and a standard error of 0.4 kJ mol−1. More information about VM6 and VM9 and the associated data sets is given in the Supporting Information, in a spreadsheet for automated generation of estimated MADs (EMADVM6 and EMADVM9) for MGCDB82. The previous extensive study with MGCDB8210 has already assessed a large number of computational chemistry methods. For our benchmark with the MG8, VM6, and VM9 models, we will thus consider just a few DFT methods that have not been assessed in that study. Where appropriate, we include dispersion corrections45 in the methods for improved overall accuracies. The methods include, on the one hand, a slightly revised form of a general-purpose method that we recently formulated (HSEB-D),15 and, on the other hand, three functionals (EDF2-D,46 WP04-D,47 and APF-D48) designed for computing specific molecular properties rather than wide ranges of thermochemical quantities. As proof-of-principle, we also use the MG8 set to formulate two simple global-hybrid DFT methods (MBh-D and PGTh) and then use VM6 and VM9 to verify the predicted accuracies. A brief overview of the DFT methods and the rationale for their inclusion in our assessment is given in the Supporting Information. The EMADs for the MG8, VM6, and VM9 models for these six methods are shown in Table 5. As mentioned above, the standard error for EMADMG8 is 0.9 kJ mol−1. We assign two standard errors, for an assumed 95% confidence interval, as the uncertainty in the EMADMG8 values. Thus, EMADMG8 values all carry an uncertainty of ±1.8 kJ mol−1. In a similar manner, the predictions with VM6 and VM9 carry uncertainties of ±2.2 and 0.8 kJ mol−1, respectively. Overall, we can see a fair, albeit not outstanding, agreement between the values estimated with the three models. With that being said, there is a general trend that tells apart the seemingly more accurate ones (e.g., HSEBD) from the ones that are less good (e.g., PGTh). In addition, in most cases, the ranges spanned by the three predicted values and their uncertainties do overlap. Thus, statistically, these models generally agree with one another, as least for the six methods examined. We take the three ranges into account and conservatively assign to each DFT method the maximum range. This yields our final EMAD values and the associated uncertainties (Table
Figure 2. Distribution of deviations for the MGCDB82 set for ωB97M-V and B97-D. E
DOI: 10.1021/acs.jctc.8b00514 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation the systems with deviations that fall within the range of ±5 kJ mol−1. This is significant as such a range roughly defines the so-called “chemical accuracy” of 1 kcal mol−1. For B97-D, the percentage for this range drops to 60%. At the other end of the spectrum, for ωB97M-V, there are ∼2% of the systems with deviations having magnitudes that are larger than 25 kJ mol−1. These are mainly associated with atomization energies in the W4-11 set52 or quantities derived from these atomization energies. In comparison, B97-D is associated with ∼10% of such large deviations. While many of these deviations are related to W4-11, there are also a substantial number of entries for other properties, notably barrier heights but also isomerization energies and noncovalent interactions. These observations show that, apart from situations in which one employs the most accurate functionals that are currently available, it is likely that the use of many DFT methods would lead to large deviations from time to time. More importantly, it may not be trivial to anticipate the circumstances under which this occurs. The ideal solution to this dilemma is perhaps the use of a highly accurate method, for example, ωB97M-V, when one performs main-group thermochemical calculations. However, such an option might not be viable, for instance, for large systems for which a screened-exchange DFT or a pure functional are the only computationally feasible choices. Thus, it has been and still is advisable to conduct benchmark studies on related model systems before proceeding with production computations on actual systems of interest. For developers of new quantum chemistry methods, the observations suggest that a highly robust method for maingroup chemistry is likely to yield an MADMGCDB82 (and, as an approximation, EMADMG8) that is smaller than ∼4 kJ mol−1. While this represents a seemingly mighty aim, we nonetheless propose it to be an initial target for the formulation and refinement of new quantum chemistry methods. A caveat of this recommendation is that, because MG8 contains a relatively small set with just 64 data points, development with this model may be prone to overfitting with a high degree of parametrization, and one should be mindful of such a potential trap. One should also keep in mind that EMADMG8 is designed to reflect the overall MAD for MGCDB82, which we select as an (presumably, realistic and robust) indicator for the accuracy for general main-group chemistry. Whether MGCDB82 is truly representative is perhaps open for debate. In this regard, we believe that the principles and methodologies proposed in the present study for devising models that contain small test sets would remain a useful tool should one wish to simplify the assessment procedure with other high-quality and extensive data sets. In the next section, we will briefly touch on this topic by examining some potential complementary tests in addition to benchmark using MG8. Beyond MG8, Prospects for a Wider Scope. The MGCDB82 set contains a diverse range of chemical properties but only for main-group molecular species, and an expanded scope would be useful. For instance, the ωB97M-V6 functional was originally formulated with MGCDB84, but in the subsequent extensive benchmark study,10 it has been tested for a prototypical transition-metal-catalyzed reaction. For that reaction, it yields fairly good accuracy, which lends further confidence to its formulation and applicability. More generally, modern DFT methods such as ωB97M-V and MN155 have increasingly shown promises for the provision of reasonable results for not just main-group molecular chemistry but also for
other properties. There is a wide range of additional chemical properties that are of fundamental and practical interests. Here we briefly consider the possible complement to MG8 with just three additional small-data-set models (and a fourth one given in the Supporting Information). The recently introduced MOR41 set51 contains 41 reactions of transition-metal complexes that are of practical relevance to typical catalytic processes. While the test set is moderately sized in terms of the number of data points, the species involved are quite large (up to ∼120 atoms), and therefore its use still involves substantial computational resources. In that study, MOR41 has been used for assessing the accuracy of almost 100 quantum chemistry methods for modeling typical homogeneous transition-metal catalysis. The raw data for the benchmark study has been made publicly available, and we will use these to attempt devising a model of smaller data sets. We note that the species in MOR41 has been intentionally limited to just largely single-reference systems, but the authors of ref 51 have hinted at an expansion to multireference species being underway. Another area in which modern general-purpose DFT methods have made good progress is the modeling of material systems. For instance, pure functionals were typically preferred over global-hybrid DFT methods in terms of the accuracy for describing materials band gaps.52 However, MN15 (a global hybrid) has shown comparable accuracy to some of the most widely used pure DFT methods in materials modeling, such as PBE,53 for the SBG31 set54 of 31 band gaps. From the perspective of computational efficiency, however, pure and screened-exchange functionals are still favored over global hybrids. In the present study, we will attempt to trim the number of data points in the SBG31 set. While this set has been used for assessing a sizable collection of DFT methods, raw data for individual deviations are at present not publicly available. We will therefore produce our own set of data for statistical analysis. We will limit ourselves to pure and screened-exchange functionals due to computational resources consideration. In addition, we will restrict the methods that we use to non-dispersion-corrected methods due to limitations in the software programs that are available to us. In this way, we have chosen 26 DFT methods (see the Supporting Information), which is by no means an exhaustive list but still reasonable when applied to a data set of 31 points. The last set that we consider is the MB08-165 set.55 It consists of 165 reactions of artificial molecules, and these species serve to “stress test” how a quantum chemistry method performs for random local minima on potential energy surfaces. The MB08-165 has been used to benchmark ∼100 methods,56 and the raw deviation data has been made publicly available. This provides a platform for statistically formulating a representative subset. We note that a smaller set of 43 reactions in MB08-165 has been included in the GMTKN55 compilation.11 In this regard, it is of interest to examine whether one can arrive at a smaller collection. The general strategy that we apply to these three sets is the same as that used on the subsets of MGCDB82. For each set, we use standard lasso and forward selection techniques to arrive at a model that contains the smallest number of data points for estimating the overall MAD as a linear combination of the individual absolute deviations. Simultaneously, we place a condition of having an R2 that is better than 0.995 for the correlation with the actual MADs for all the assessed quantum chemistry methods. In this manner, the numbers of data points F
DOI: 10.1021/acs.jctc.8b00514 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
complete substitutes for their parent sets. Small models are useful for rapid initial screening of theoretical models, but they do not provide details at the same level as those offered by data sets that are more comprehensive. In this regard, the richness of information contained in large data sets remains an invaluable tool for computational chemists.
in MOR41, SBG31, and MB08-165 are reduced to 13, 5, and 13, respectively. For our discussion, we will refer to these models of trimmed sets as MOR13, SBG5, and MB13. The reactions contained in MOR13 are abbreviated and shown in Figure 3. They are numbered according to ref 51. Within the
■
CONCLUDING REMARKS In the present study, we have employed literature data to survey the deviations for a wide range of computational chemistry methods for large test sets. These include mainly a major collection of subsets of MGCDB84 (that we call MGCDB82), while we have also briefly examined the GMTKN55 compendia (see the Supporting Information). Our analysis shows that the overall statistical indicators for MGCDB82 and GMTKN55 correlate reasonably well with one another. We also find different statistical metrics such as mean absolute deviations (MADs) and mean absolute fractional deviations provide similar qualitative observations in terms of categorizing theoretical methods for their overall accuracies. We use the MGCDB82 set to devise models of collections of smaller sets for rapid evaluation of new quantum chemistry procedures. An aim is for the model to provide good estimated MADs (EMADs) for the overall MGCDB82 set. Another objective is to yield qualitative insights into the performance for the eight general categories of chemical properties contained in MGCDB82. The approaches that we have used involve standard statistical model-building tools such as lasso regularization and forward selection techniques, in addition to comparison based on straightforward linear regressions. As a result, we have devised the MG8 model. It contains a total of just 64 data points spreading across eight downsized subsets, each representing a category in MGCDB82. The EMADs obtained with MG8 (EMADMG8) are a weighted linear combination of the MADs for the underlying subsets. In a similar manner, each of these eight MADs is estimated by linearly combining the data points (absolute deviations) in the set. The EMADMG8 values correlate with the actual MADs for MGCDB82 for all the (∼200) methods assessed in ref 10 with an R2 of 0.988. While the use of EMADMG8 represents an economical means for assessing the general accuracy of quantum chemistry methods, we caution that additional testing should not be completely foregone. This is because, for the large number of DFT methods previously assessed with MGCDB82, an “average” functional is often associated with a significant number of outliers. Importantly, it is not straightforward to foretell the systems for which such methods do not perform well. With this caveat in mind, we propose a EMADMG8 value of ∼4 kJ mol−1 to be a sensible (though challenging) target for quantum chemistry method development. In relation to the matter of diversity in chemical properties, it is desirable to widen the scope of MG8 (as an approximation to the full MBCDB82 set), which contains only main-group molecular species. Thus, we have formulated additional smalldata-set models for a higher degree of diversity. They include the MOR13 model for metal−organic reactions, the SBG5 model for semiconductor band gaps, and MB13 for stresstesting methods with artificial species. Assisting tools for carrying out necessary computations for the range of models and for analyzing the results are given in the Supporting Information.
Figure 3. Reactions contained in the MOR13 model. Ligands in the transition-metal complexes that do not participate directly are omitted for clarity. The reactions are numbered according to ref 51. The abbreviations for chemical groups and species are “Me” for methyl, “Ac” for acetyl”, “Cy” for cyclohexyl, “dmpe” for “1,2-bis(dimethylphosphino)ethane”, and “Ph” for phenyl.
MOR41 set, the numbers of reactions that involve 3d, 4d, and 5d transition metals are 10, 21, and 10, respectively. The corresponding numbers in MOR13 are 3, 8, and 2. For SBG5, the systems that it contains are AlSb, BaTe, C (diamond), CdSe, and ZnS. Let us use HSEB-D as a case study for the use of MOR13, SBG5, and MB13 for assessing a method trained with maingroup molecular properties (the E3 set24−26 in this case). The performance of HSEB(-D) for SBG31 and MB08-165 has been previously examined, and the MADs for these two sets are 0.33 eV and 25.0 kJ mol−1, respectively. In comparison, the corresponding EMADSBG5 and EMADMB13 values are 0.33 eV and 26.3 kJ mol−1. The agreements provide another piece of evidence supporting the use of models of statistically generated small test sets. What about transition-metal chemistry for which the accuracy of HSEB-D has yet to be examined? With the MOR13 model, the EMAD is 23.2 kJ mol−1. This is far from being in the pack of the most accurate methods examined in ref 67 (e.g., ωB97X-V57 has an MAD of 9.2 kJ mol−1). However, it is among the methods with above-average accuracies, which is the conclusion that we have come to on the basis of our original analysis15 as well as those performed in the present study with MG8. At this point, a relevant question is whether to incorporate MOR13, SBG5, and MB13 into the MG8 model. We believe in the utility of a single well-balanced metric for the overall accuracy for a diverse range of systems and properties. However, we do not consider a conclusion has been reached regarding the definition of “well-balanced”. We therefore perceive that, at this stage, it is perhaps premature for us to integrate these additional tests, and we recommend their use simply as independent assessments. We also reiterate that we do not advocate the use of MG8 and other small models as G
DOI: 10.1021/acs.jctc.8b00514 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
■
(10) Mardirossian, N.; Head-Gordon, M. Thirty Years of Density Functional Theory in Computational Chemistry: An Overview and Extensive Assessment of 200 Density Functionals. Mol. Phys. 2017, 115, 2315−2372. (11) Goerigk, L.; Hansen, A.; Bauer, C. A.; Ehrlich, S.; Najibi, A.; Grimme, S. A Look at the Density Functional Theory Zoo with the Advanced GMTKN55 Database for General Main Group Thermochemistry, Kinetics and Noncovalent Interactions. Phys. Chem. Chem. Phys. 2017, 19, 32184−32215. (12) Chan, B.; Gilbert, A. T. B.; Gill, P. M. W.; Radom, L. Performance of Density Functional Theory Procedures for the Calculation of Proton-Exchange Barriers: Unusual Behavior of M06Type Functionals. J. Chem. Theory Comput. 2014, 10, 3777−3783. (13) Chan, B.; Song, J. W.; Kawashima, Y.; Hirao, K. Toward the Complete Range Separation of Non-Hybrid Exchange−Correlation Functional. J. Comput. Chem. 2015, 36, 871−877. (14) Chan, B.; Song, J. W.; Kawashima, Y.; Hirao, K. Performance of the OP Correlation Functional in Relation to its Formulation: Influence of the Exchange Component and the Effect of Incorporating Same-Spin Correlations. J. Comput. Chem. 2016, 37, 1306−1312. (15) Chan, B.; Kawashima, Y.; Hirao, K. Correlation Functional in Screened-Exchange Density Functional Theory Procedures. J. Comput. Chem. 2017, 38, 2307−2315. (16) Chan, B.; Kawashima, Y.; Hirao, K. The reHISS Three-Range Exchange Functional with an Optimal Variation of Hartree−Fock and Its Use in the reHISSB-D Density Functional Theory Method. J. Comput. Chem. DOI: 10.1002/jcc.25383. (17) Chan, B.; Coote, M. L.; Radom, L. G4-SP, G4(MP2)-SP, G4sc, and G4(MP2)-sc: Modifications to G4 and G4(MP2) for the Treatment of Medium-Sized Radicals. J. Chem. Theory Comput. 2010, 6, 2647−2653. (18) Chan, B.; Deng, J.; Radom, L. G4(MP2)-6X: A Cost-Effective Improvement to G4(MP2). J. Chem. Theory Comput. 2011, 7, 112− 120. (19) Chan, B.; Radom, L. W1X-1 and W1X-2: W1-Quality Accuracy with an Order of Magnitude Reduction in Computational Cost. J. Chem. Theory Comput. 2012, 8, 4259−4269. (20) Chan, B.; Radom, L. W3X: A Cost-Effective Post-CCSD(T) Composite Procedure. J. Chem. Theory Comput. 2013, 9, 4769−4778. (21) Chan, B.; Radom, L. W2X and W3X-L: Cost-Effective Approximations to W2 and W4 with kJ mol−1 Accuracy. J. Chem. Theory Comput. 2015, 11, 2109−2119. (22) Chan, B. Unification of the W1X and G4(MP2)-6X Composite Protocols. J. Chem. Theory Comput. 2017, 13, 2642−2649. (23) Chan, B. How to Computationally Calculate Thermochemical Properties Objectively, Accurately, and as Economically as Possible. Pure Appl. Chem. 2017, 89, 699−713. (24) Chan, B.; Radom, L. Obtaining Good Performance with Tripleζ-Type Basis Sets in Double-Hybrid Density Functional Theory Procedures. J. Chem. Theory Comput. 2011, 7, 2852−2863. (25) Chan, B.; Gill, P. M. W.; Radom, L. Performance of GradientCorrected and Hybrid Density Functional Theory: Role of the Underlying Local Density Approximation and the Gradient Correction. J. Chem. Theory Comput. 2012, 8, 4899−4906. (26) Chan, B.; Karton, A.; Raghavachari, K.; Radom, L. RestrictedOpen-Shell G4(MP2)-Type Procedures. J. Phys. Chem. A 2016, 120, 9299−9304. (27) Lynch, B. J.; Truhlar, D. G. Small Representative Benchmarks for Thermochemical Calculations. J. Phys. Chem. A 2003, 107, 8996− 8999. (28) Zheng, J.; Zhao, Y.; Truhlar, D. G. The DBH24/08 Database and Its Use to Assess Electronic Structure Model Chemistries for Chemical Reaction Barrier Heights. J. Chem. Theory Comput. 2009, 5, 808−821. (29) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Petersson, G. A.; Nakatsuji, H.; Li, X.; Caricato, M.; Marenich, A. V.; Bloino, J.; Janesko, B. G.; Gomperts, R.; Mennucci, B.; Hratchian, H. P.; Ortiz, J. V.; Izmaylov, A. F.; Sonnenberg, J. L.; Williams-Young, D.; Ding, F.;
ASSOCIATED CONTENT
S Supporting Information *
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.8b00514.
■
Additional discussions (PDF) Templates for obtaining EMADs with these models (XLSX) Geometries for the species in the data sets of MG8, VM6, VM9, MOR13, SBG5, and MB13 discussed in the main text and those for XE6 discussed in the Supporting Information (ZIP)
AUTHOR INFORMATION
Corresponding Author
*E-mail:
[email protected]. ORCID
Bun Chan: 0000-0002-0082-5497 Funding
We gratefully acknowledge research funding from the Japan Society for the Promotion of Science (JSPS) (Grant number 16H07074001 and Program for Advancing Strategic International Networks to Accelerate the Circulation of Talented Researchers) and generous grants of computer time from the RIKEN Advanced Center for Computing and Communication (ACCC) (Project Q18266), Japan, and National Computational Infrastructure (NCI), Australia. Notes
The author declares no competing financial interest.
■
REFERENCES
(1) Sholl, D.; Steckel, J. A. Density Functional Theory: A Practical Introduction; Wiley: Hoboken, 2009; DOI: 10.1002/9780470447710. (2) Boese, A. D.; Martin, J. M. L. Development of Density Functionals for Thermochemical Kinetics. J. Chem. Phys. 2004, 121, 3405−3416. (3) Zhao, Y.; Truhlar, D. G. The M06 Suite of Density Functionals for Main Group Thermochemistry, Thermochemical Kinetics, Noncovalent Interactions, Excited States, and Transition Elements: Two New Functionals and Systematic Testing of Four M06-Class Functionals and 12 Other Functionals. Theor. Chem. Acc. 2008, 120, 215−241. (4) Chai, J.-D.; Head-Gordon, M. Long-Range Corrected Hybrid Density Functionals with Damped Atom-Atom Dispersion Corrections. Phys. Chem. Chem. Phys. 2008, 10, 6615−6620. (5) Yu, H. S.; He, X.; Li, S. L.; Truhlar, D. G. MN15: A Kohn-Sham Global-Hybrid Exchange-Correlation Density Functional with Broad Accuracy for Multi-Reference and Single-Reference Systems and Noncovalent Interactions. Chem. Sci. 2016, 7, 5032−5051. (6) Mardirossian, N.; Head-Gordon, M. ωB97M-V: A Combinatorially Optimized, Range-Separated Hybrid, Meta-GGA Density Functional with VV10 Nonlocal Correlation. J. Chem. Phys. 2016, 144, 214110−1−23. (7) Curtiss, L. A.; Raghavachari, K.; Redfern, P. C.; Pople, J. A. Assessment of Gaussian-2 and Density Functional Methods for the Computation of Enthalpies of Formation. J. Chem. Phys. 1997, 106, 1063−1079. (8) Curtiss, L. A.; Redfern, P. C.; Raghavachari, K.; Pople, J. A. Assessment of Gaussian-2 and Density Functional Methods for the Computation of Ionization Energies and Electron Affinities. J. Chem. Phys. 1998, 109, 42−55. (9) Becke, A. D. Density-Functional Thermochemistry. V. Systematic Optimization of Exchange−Correlation Functionals. J. Chem. Phys. 1997, 107, 8554−8560. H
DOI: 10.1021/acs.jctc.8b00514 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
(38) Sall, J.; Lehman, A.; Stephens, M.; Loring, S. JMP Start Statistics: A Guide to Statistics and Data Analysis Using JMP, Sixth ed.; SAS Institute Inc.: Cary NC, 2017. (39) Ŕ ezác,̌ J.; Riley, K. E.; Hobza, P. S66: A Well-balanced Database of Benchmark Interaction Energies Relevant to Biomolecular Structures. J. Chem. Theory Comput. 2011, 7, 2427−2438. (40) Lao, K. U.; Herbert, J. M. Accurate and Efficient Quantum Chemistry Calculations for Noncovalent Interactions in Many-Body Systems: The XSAPT Family of Methods. J. Phys. Chem. A 2015, 119, 235−252. (41) Bauzá, A.; Alkorta, I.; Frontera, A.; Elguero, J. On the Reliability of Pure and Hybrid DFT Methods for the Evaluation of Halogen, Chalcogen, and Pnicogen Bonds Involving Anionic and Neutral Electron Donors. J. Chem. Theory Comput. 2013, 9, 5201− 5210. (42) Yuan, Y.; Mills, M. J. L.; Popelier, P. L. A.; Jensen, F. Comprehensive Analysis of Energy Minima of the 20 Natural Amino Acids. J. Phys. Chem. A 2014, 118, 7876−7891. (43) Manna, D.; Martin, J. M. L. What Are the Ground State Structures of C20 and C24? An Explicitly Correlated Ab Initio Approach. J. Phys. Chem. A 2016, 120, 153−160. (44) Karton, A.; Daon, S.; Martin, J. M. L. W4−11: A HighConfidence Benchmark Dataset for Computational Thermochemistry Derived from First-Principles W4 Data. Chem. Phys. Lett. 2011, 510, 165−178. (45) Grimme, S. Density Functional Theory with London Dispersion Corrections. WIREs Comput. Mol. Sci. 2011, 1, 211−228. (46) Lin, C. Y.; George, M. W.; Gill, P. M. W. EDF2: A Density Functional for Predicting Molecular Vibrational Frequencies. Aust. J. Chem. 2004, 57, 365−370. (47) Wiitala, K. W.; Hoye, T. R.; Cramer, C. J. Hybrid Density Functional Methods Empirically Optimized for the Computation of 13 C and 1H Chemical Shifts in Chloroform Solution. J. Chem. Theory Comput. 2006, 2, 1085−1092. (48) Austin, A.; Petersson, G.; Frisch, M. J.; Dobek, F. J.; Scalmani, G.; Throssell, K. A Density Functional with Spherical Atom Dispersion Terms. J. Chem. Theory Comput. 2012, 8, 4989−5007. (49) Henderson, T. M.; Janesko, B. G.; Scuseria, G. E. Generalized Gradient Approximation Model Exchange Holes for Range-Separated Hybrids. J. Chem. Phys. 2008, 128, 194105−1−9. (50) Grimme, S. Semiempirical GGA-Type Density Functional Constructed with a Long-Range Dispersion Correction. J. Comput. Chem. 2006, 27, 1787−1799. (51) Dohm, S.; Hansen, A.; Steinmetz, M.; Grimme, S.; Checinski, M. P. Comprehensive Thermochemical Benchmark Set of Realistic Closed-Shell Metal Organic Reactions. J. Chem. Theory Comput. 2018, 14, 2596−2608. (52) Janesko, B. G.; Henderson, T. M.; Scuseria, G. E. Screened Hybrid Density Functionals for Solid-State Chemistry and Physics. Phys. Chem. Chem. Phys. 2009, 11, 443−454. (53) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865−3868. (54) Peverati, R.; Truhlar, D. G. Performance of the M11-L Density Functional for Bandgaps and Lattice Constants of Unary and Binary Semiconductors. J. Chem. Phys. 2012, 136, 134704−1−10. (55) Korth, M.; Grimme, S. Mindless” DFT Benchmarking. J. Chem. Theory Comput. 2009, 5, 993−1003. (56) Goerigk, L.; Grimme, S. A Thorough Benchmark of Density Functional Methods for General Main Group Thermochemistry, Kinetics, And Noncovalent Interactions. Phys. Chem. Chem. Phys. 2011, 13, 6670−6688. (57) Mardirossian, N.; Head-Gordon, M. ωB97X-V: A 10Parameter, Range-Separated Hybrid, Generalized Gradient Approximation Density Functional with Nonlocal Correlation, Designed by a Survival-of-the-Fittest Strategy. Phys. Chem. Chem. Phys. 2014, 16, 9904−9924.
Lipparini, F.; Egidi, F.; Goings, J.; Peng, B.; Petrone, A.; Henderson, T.; Ranasinghe, D.; Zakrzewski, V. G.; Gao, J.; Rega, N.; Zheng, G.; Liang, W.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Throssell, K.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M. J.; Heyd, J. J.; Brothers, E. N.; Kudin, K. N.; Staroverov, V. N.; Keith, T. A.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A. P.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Millam, J. M.; Klene, M.; Adamo, C.; Cammi, R.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Farkas, O.; Foresman, J. B.; Fox, D. J. Gaussian 16, Revision A.03; Gaussian, Inc.: Wallingford, CT, 2016. (30) Shao, Y.; Gan, Z.; Epifanovsky, E.; Gilbert, A. T. B.; Wormit, M.; Kussmann, J.; Lange, A. W.; Behn, A.; Deng, J.; Feng, X.; Ghosh, D.; Goldey, M.; Horn, P. R.; Jacobson, L. D.; Kaliman, I.; Khaliullin, R. Z.; Kús, T.; Landau, A.; Liu, J.; Proynov, E. I.; Rhee, Y. M.; Richard, R. M.; Rohrdanz, M. A.; Steele, R. P.; Sundstrom, E. J.; Woodcock, H. L., III; Zimmerman, P. M.; Zuev, D.; Albrecht, B.; Alguire, E.; Austin, B.; Beran, G. J. O.; Bernard, Y. A.; Berquist, E.; Brandhorst, K.; Bravaya, K. B.; Brown, S. T.; Casanova, D.; Chang, C.M.; Chen, Y.; Chien, S. H.; Closser, K. D.; Crittenden, D. L.; Diedenhofen, M.; DiStasio, R. A., Jr.; Dop, H.; Dutoi, A. D.; Edgar, R. G.; Fatehi, S.; Fusti-Molnar, L.; Ghysels, A.; Golubeva-Zadorozhnaya, A.; Gomes, J.; Hanson-Heine, M. W. D.; Harbach, P. H. P.; Hauser, A. W.; Hohenstein, E. G.; Holden, Z. C.; Jagau, T.-C.; Ji, H.; Kaduk, B.; Khistyaev, K.; Kim, J.; Kim, J.; King, R. A.; Klunzinger, P.; Kosenkov, D.; Kowalczyk, T.; Krauter, C. M.; Lao, K. U.; Laurent, A.; Lawler, K. V.; Levchenko, S. V.; Lin, C. Y.; Liu, F.; Livshits, E.; Lochan, R. C.; Luenser, A.; Manohar, P.; Manzer, S. F.; Mao, S.-P.; Mardirossian, N.; Marenich, A. V.; Maurer, S. A.; Mayhall, N. J.; Oana, C. M.; Olivares-Amaya, R.; O’Neill, D. P.; Parkhill, J. A.; Perrine, T. M.; Peverati, R.; Pieniazek, P. A.; Prociuk, A.; Rehn, D. R.; Rosta, E.; Russ, N. J.; Sergueev, N.; Sharada, S. M.; Sharmaa, S.; Small, D. W.; Sodt, A.; Stein, T.; Stück, D.; Su, Y.-C.; Thom, A. J. W.; Tsuchimochi, T.; Vogt, L.; Vydrov, O.; Wang, T.; Watson, M. A.; Wenzel, J.; White, A.; Williams, C. F.; Vanovschi, V.; Yeganeh, S.; Yost, S. R.; You, Z.-Q.; Zhang, I. Y.; Zhang, X.; Zhou, Y.; Brooks, B. R.; Chan, G. K. L.; Chipman, D. M.; Cramer, C. J.; Goddard, W. A., III; Gordon, M. S.; Hehre, W. J.; Klamt, A.; Schaefer, H. F., III; Schmidt, M. W.; Sherrill, C. D.; Truhlar, D. G.; Warshel, A.; Xua, X.; Aspuru-Guzik, A.; Baer, R.; Bell, A. T.; Besley, N. A.; Chai, J.-D.; Dreuw, A.; Dunietz, B. D.; Furlani, T. R.; Gwaltney, S. R.; Hsu, C.-P.; Jung, Y.; Kong, J.; Lambrecht, D. S.; Liang, W.; Ochsenfeld, C.; Rassolov, V. A.; Slipchenko, L. V.; Subotnik, J. E.; Van Voorhis, T.; Herbert, J. M.; Krylov, A. I.; Gill, P. M. W.; Head-Gordon, M. Advances in Molecular Quantum Chemistry Contained in the QChem 4 Program Package. Mol. Phys. 2015, 113, 184−215. (31) (a) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H−Pu. J. Chem. Phys. 2010, 132, 154104. (b) See also: https://www.chemie. uni-bonn.de/pctc/mulliken-center/software/dft-d3/ (accessed Jan 2018). (32) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. A. J. Chem. Phys. 1980, 72, 650−654. (33) Weigend, F.; Ahlrichs, R. Balanced Basis Sets of Split Valence, Triple Zeta Valence and Quadruple Zeta Valence Quality for H to Rn: Design and Assessment of Accuracy. Phys. Chem. Chem. Phys. 2005, 7, 3297−3305. (34) Heyd, J.; Peralta, J. E.; Scuseria, G. E.; Martin, R. L. Energy Band Gaps and Lattice Parameters Evaluated with the Heyd-ScuseriaErnzerhof Screened Hybrid Functional. J. Chem. Phys. 2005, 123, 174101. (35) Jones, B.; Sall, J. JMP Statistical Discovery Software. WIREs Comp. Stat. 2011, 3, 188−194. (36) JMP, Version 14; SAS Institute Inc.: Cary NC, 2018. (37) Mardirossian, N.; Head-Gordon, M. Survival of the Most Transferable at the Top of Jacob’s Ladder: Defining and Testing the ωB97M(2) Double Hybrid Density Functional. J. Chem. Phys. 2018, 148, 241736. I
DOI: 10.1021/acs.jctc.8b00514 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX