J. Chem. Inf. Comput. Sci. 1998, 38, 639-645
639
Prediction of Critical Temperatures and Pressures of Industrially Important Organic Compounds from Molecular Structure Brian E. Turner, Chandra L. Costello, and Peter C. Jurs* Department of Chemistry, Penn State University, 152 Davey Laboratory, University Park, Pennsylvania 16802 Received January 23, 1998
Quantitative-structure property relationships methods are used to develop mathematical models to predict critical temperatures and pressures of a diverse set of organic compounds taken from the Design Institute for Physical Property Data (DIPPR) database. Each compound is represented with calculated molecular structure descriptors that encode its topological, electronic, geometrical, and other features. Subsets of descriptors are selected with simulated annealing and genetic algorithms. Models to predict the critical properties are constructed using multiple linear regression analysis and computational neural networks with errors comparable to the experimental errors of the critical property data. INTRODUCTION
Quantitative structure-property relationship (QSPR) methods have been used to develop models for the prediction of two fundamental physical properties of organic compounds: critical temperature and pressure. The critical temperature is the temperature above which a gas cannot be liquefied, and the minimum pressure required for liquefaction at this temperature is the critical pressure. Knowledge of these critical properties for organic compounds is very important to industrial chemical engineers. The values of critical temperatures and pressures are needed in equation of state calculations of thermodynamic and transport properties. In addition, these physical properties are important in high-pressure phase equilibrium processes. Such data are needed, for example, for studies of enhanced oil recovery and supercritical fluid extraction.1 As a result of the importance of these properties, a number of different approaches have been taken in order to calculate the values of critical temperatures and pressures. One common approach is to use knowledge of other physical properties. An example of this method is the estimation of the critical point using the Peng-Robinson equation of state.1 This method uses the van der Waals volume, the normal boiling point, and the acentric factor of a compound as variables. A drawback of such a method is the reliance on experimental data which can be expensive or time-consuming to gather. An alternative approach is to build models to predict critical temperatures of organic compounds using structural information alone.2 This method employs quantitative structure-property relationships (QSPRs) to develop its predictions. Different methodologies have been reported to develop QSPRs. One such approach has been through the use of the Automated Data Analysis and Pattern recognition Toolkit (ADAPT).3 This interactive software system utilizes multiple linear regression analysis (MLR) and computational neural networks (CNN) as model-building routines to produce QSPRs. A goal of the present study is to use these
methods to improve on the prediction capabilities using QSPR for critical temperatures. In general, development of a QSPR involves three steps: structural encoding, feature selection, and model building. Structural encoding involves the use of numerical descriptors to encode the structural features of a compound. Feature selection is then employed to determine which subset of the descriptors best relates to the property of interest. Models built from the best subset of descriptors form a direct link between descriptors and the property of interest. Finally, validation determines the level of the model’s predictive capabilities for unknown compounds. Recently, several changes in the ADAPT software have lead to the belief that improvements have been made in the QSPR methodology.4-6 Several new approaches for numerically encoding a compound’s structure have been implemented. For example, development of carbon type descriptors allows for structural descriptors such as the number of secondary sp3 carbons to be determined, and hydrogenbonding descriptors encode the ability of compounds to engage in hydrogen bonding interactions. The development of better descriptors is a continuing goal of QSPR studies. In addition to new descriptors, the ADAPT software system has had updates in its approach to subjective feature selection.7,8 Previously, the method of leaps-and-bounds regression was applied in order to build regression models. However, new developments have allowed a genetic algorithm and simulated annealing to replace the old feature selection methods. These new means of feature selection improve results in QSPR models because of their capability to systematically search through a larger number of descriptors and discard from consideration any variables that do not contribute to the model being generated. The development of new models for predicting critical temperature is the first step here. Once assessments of the methodology are made, the second goal of the study is to develop a predictive equation for critical pressures using the same set of compounds. Together, the two equations can support the prediction of the critical point for various types of organic compounds.
S0095-2338(98)00005-5 CCC: $15.00 © 1998 American Chemical Society Published on Web 05/20/1998
640 J. Chem. Inf. Comput. Sci., Vol. 38, No. 4, 1998 METHODOLOGY
The Data Set. Critical temperatures and pressures were obtained for a set of 165 compounds from the DIPPR database.9 The functional groups found among the data set compounds include alcohols, ketones, esters, carboxylic acids, aldehydes, phenols, ethers, nitriles, and amines. The hydrocarbons include alkanes, alkenes, benzenes, cycloalkanes, and cycloalkenes. A complete list of the compounds in the data set and their corresponding experimental critical temperatures and pressures is shown in Table 1. The data set contained 76 hydrocarbons, 61 oxygen-containing compounds, 23 nitrogen-containing compounds, and 5 halogencontaining compounds. The compounds range in size from propylene (molecular weight of 42) to quinoline (molecular weight of 129). For the linear regression studies, the data set was divided into a training set (tset) of 147 compounds for model building and an external prediction set (pset) of 18 compounds for model validation. Both subsets were chosen to ensure that a diverse set of compounds was present. For the computational neural network (CNN) studies, a crossvalidation set (cvset) of 15 compounds was chosen, leaving 132 compounds in the tset, and the pset remained the same. The cvset is a subset of compounds used to help find an optimal set of weights and biases during CNN training, and it is also used to avoid overtraining of the CNN. The range of values for both critical temperatures and pressures was large, with the critical temperatures ranging from 365 to 782 K and the critical pressures ranging from 11.94 to 55.47 atm. The experimental uncertainties associated with the values used in this data set are determined by the members of the DIPPR staff. Calculations show that the experimental critical temperatures and pressures in this study have a 1.3% and 1.5% error, respectively. These uncertainties indicate that the ideal errors associated with our models should be approximately 8 K and 0.5 atm. Molecular Modeling. The 165 compounds were sketched into HyperChem10 and transferred to the workstation in the form of connection tables. A more extensive energy minimization was done using MOPAC11 with the PM3 Hamiltonian.12 The resulting three-dimensional conformations were necessary so that they could be used later to calculate geometrical descriptors. Descriptor Generation. Once the structures were entered, molecular structure descriptors were generated. Three categories of descriptors were calculated: topological, geometrical, and electronic. Topological descriptors include values such as the total number of atoms, the number of double bonds, molecular connectivity indices, and path counts. Geometrical descriptors encode structural features such as the surface area of a compound, its moments of inertia, and its molecular volume. Electronic descriptors provide information about such features as the charge on various atoms within the molecule or dipole moment. In addition to these descriptors, hybrids known as charged partial surface area (CPSA) descriptors use both geometrical and electronic information to encode the compounds’ ability to engage in polar interactions.13 Finally, other hybrid descriptors account for the effects of hydrogen bonding. A total of 192 descriptors was calculated for each compound in the data set.
TURNER
ET AL.
Descriptor Screening. The process of descriptor generation produced a large number of descriptors for each compound. However, to reduce the likelihood of chance correlations being found during model development, the number of descriptors had to be reduced to be less than 60% of the number of observations. The reduction of the number of descriptors involved a series of steps. First, all descriptors with 90% identical values were removed because they did not offer valuable information that allowed for discrimination between compounds. Next, when pairs of descriptors were found to have pairwise correlations greater than 0.95, one of them was removed, while the other remained in order to maintain a small pool of information rich descriptors. Finally, the descriptor pool was reduced using vector space descriptor analysis (VSDA)14 by ranking the descriptors based on an orthogonalization procedure. The descriptors with the lowest ranking were removed until the pool was small enough. Regression Analysis. After a reduced pool of descriptors was developed, multiple linear regression analysis was used to develop a model that predicted the property of interest:
Y ) b0 + b1X1 + b2X2 + ... + bnXn In this equation, Y is the value of the physical property of a compound, b0 is the y-intercept of the regression model, b1 through bn are the various coefficients for the descriptors determined by the regression model, and X1 through Xn are the specific descriptor values for the compound. To formulate regression models, two methods of feature selection were utilized. A genetic algorithm and simulated annealing allowed for regression models of various sizes to be built from the larger pool of descriptors. Models of different sizes were investigated until the inclusion of an additional descriptor did not significantly improve the model statistics. Once developed, an attempt to improve the model was done using an interactive regression analysis (IRA) routine. This program allowed individual, specific descriptors to be added or removed in order to seek the optimum model statistics. When adding or changing a descriptor no longer improved the rms error of the model, it was determined that the best subset had been found. Once a model was developed, regression diagnostics were performed in order to determine if any compounds were statistical outliers. Examples of these tests include studentized residual, leverage value, and Cook’s distance.15 Four out of six of these statistical tests must fail in order for a compound to be classified as an outlier. If an outlier is found, the compound is dropped from the tset, and the model coefficients are recalculated. Model Validation. Once a regression model was selected, several investigations were conducted to determine if the model was statistically valid. The first thing examined was basic model statistics such as F, T, and p-values. Next, the model was tested for multicollinearities using the variance inflation factor (VIF).16 A VIF greater than 10 for a descriptor generally indicates the presence of an unacceptably large multicollinearity with other descriptors. Graphical approaches were then used to search for bias. To do this, experimental versus predicted plots and residual plots were constructed and examined for trends. Finally, an external pset, which was not used during feature selection or
PREDICTION
OF
CRITICAL TEMPERATURES
AND
PRESSURES
J. Chem. Inf. Comput. Sci., Vol. 38, No. 4, 1998 641
Table 1. Experimental and Predicted Values of Each Compound for All Four Models 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73
compound name
exptl pressure
linear predicted
CNN model predicted
exptl temp
linear predicted
CNN model predicted
3-chloropropene propylene allyl alcoholb ethyl formate methyl acetate propionic acida n-propyl chloride propaneb methyl ethyl ether propanol methylal isopropylamine n-propylamine trimethylamine 1,3-butadiene n-butyronitrile 1-butene cis-2-butene trans-2-butene isobutene methyl ethyl ketone tetrahydrofuran n-butyric acid ethyl acetate isobutyric acid methyl propionate n-propyl formateb n-butane isobutane n-butanol sec-butyl alcohol tert-butyl alcohol diethyl ether isobutyl alcohol n-butylamine diethylamine pyridine cyclopentane 1-pentene diethyl ketone 2-pentanone ethyl propionate isobutyl formate n-propyl acetate valeric acid isopentane neopentane n-pentane 2-methyl-2-butanol 3-methyl-1-butanol 1-pentanol bromobenzene chlorobenzene benzene phenola,c aniline 2-methylpyridine cyclohexaneb cyclohexanone cyclohexane 1-hexene methylcyclopentane cyclohexanol 2-hexanone methyl isobutyl ketone ethyl n-butyrate isobutyl acetateb 2,2-dimethylbutane 2,3-dimethylbutane n-hexaneb 2-methylpentane 3-methylpentane diisopropyl ether
46.48 45.52 55.47 46.80 46.88 45.57 45.20 41.92 43.42 51.02 39.00 44.80 46.78 40.20 42.73 37.40 39.67 41.51 40.49 39.47 41.00 51.20 40.11 38.29 40.00 39.52 40.10 37.47 36.00 43.65 41.39 39.20 35.92 42.39 41.45 36.60 55.60 44.43 34.83 36.90 36.46 33.18 38.30 33.16 37.60 33.37 31.57 33.26 38.29 38.29 38.29 44.60 44.60 48.34 60.50 52.40 43.23 42.93 38.00 40.22 30.99 37.35 37.00 32.80 32.30 29.11 29.71 30.40 30.86 29.85 29.71 30.83 28.42
46.66 44.40 50.55 43.94 44.75 48.28 43.66 40.26 42.51 47.53 41.34 47.24 44.31 45.43 44.43 43.32 39.96 39.67 38.89 43.71 41.58 47.91 42.24 38.05 43.25 39.40 40.24 36.20 37.18 43.75 43.82 39.07 37.27 42.76 40.32 37.53 51.18 41.99 36.00 37.03 36.27 34.14 34.93 34.77 36.91 33.27 29.24 32.71 37.74 39.55 40.07 45.54 45.20 48.27 NA 45.54 43.57 39.94 40.26 37.41 32.37 35.89 41.55 31.98 32.35 28.04 29.06 29.34 31.71 29.48 28.98 30.47 30.57
46.72 45.68 52.38 43.66 44.27 47.56 44.33 40.96 42.15 48.95 41.12 45.32 45.93 39.39 44.95 41.85 40.77 39.93 39.85 41.41 40.90 50.15 41.73 37.49 41.29 39.88 39.67 36.96 35.12 44.58 41.85 39.41 36.45 41.43 41.84 38.88 52.08 45.14 35.89 37.03 35.52 34.02 34.26 33.99 36.48 32.91 30.60 33.48 37.49 37.83 39.88 46.39 45.46 49.44 NA 48.33 43.50 41.51 39.38 38.60 31.36 36.35 37.12 31.27 32.11 29.57 29.67 31.02 31.00 30.29 29.81 31.17 30.85
514.5 364.8 545.0 508.4 506.8 604.0 503.1 369.8 437.8 536.7 480.6 471.8 496.9 433.3 425.4 582.3 419.6 435.6 428.6 417.9 535.5 540.1 628.0 523.3 609.1 530.6 538.0 425.2 408.1 562.9 536.0 506.2 466.7 547.7 531.9 496.6 620.0 511.8 464.8 561.0 561.1 546.0 551.3 549.4 651.0 460.4 433.8 469.6 545.1 579.5 586.1 670.1 632.3 562.2 694.3 699.0 621.0 560.4 629.1 553.5 504.0 532.8 625.1 587.0 571.4 571.0 561.0 488.8 500.0 507.4 497.5 504.4 500.0
505.8 381.5 532.4 514.5 492.7 600.1 495.7 366.5 432.4 521.6 487.3 465.1 494.3 443.3 449.1 595.4 433.4 430.4 430.2 417.0 522.1 520.9 626.6 514.2 621.7 531.3 543.5 423.9 390.3 557.6 546.9 509.9 460.4 534.3 531.8 490.4 614.6 519.4 471.8 562.2 549.0 549.1 561.3 561.4 648.8 459.6 432.4 464.7 548.3 579.5 591.5 679.8 649.7 560.0 681.2 679.5 635.8 554.3 635.8 551.6 505.7 547.8 635.7 580.4 567.0 568.6 548.8 494.1 497.0 500.3 490.1 504.9 507.4
508.8 389.5 538.5 516.1 503.6 592.0 501.5 384.6 434.2 531.5 491.2 467.1 501.7 438.1 436.2 595.2 427.7 424.3 424.7 410.4 523.2 525.0 623.3 527.8 612.8 532.8 538.5 421.7 398.0 554.5 550.4 506.8 460.9 541.1 529.6 494.1 607.9 515.7 467.5 556.2 548.7 551.0 552.3 552.3 649.7 457.5 428.4 462.3 548.8 573.0 581.6 672.8 643.6 558.1 689.4 688.2 635.1 554.5 632.5 553.3 507.4 547.1 632.0 576.8 564.2 572.1 552.6 493.6 499.0 503.9 491.1 506.0 508.0
642 J. Chem. Inf. Comput. Sci., Vol. 38, No. 4, 1998
TURNER
ET AL.
Table 1 (Continued) compound name 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146
di-n-propyl ether diisopropylamine di-n-propylamine benzaldehyde toluene m-cresol o-cresol m-toluidine o-toluidine ethylcyclopentane 2,3-dimethylpenane n-heptane 2-methylhexane 3-methylhexane ethylbenzene m-xylene p-xylene 2,6-xylenol N,N-dimethylaniline cis-1,2-dimethylcyclohexane trans-1,2-dimethylcyclohexane cis-1,3-dimethylcyclohexane trans-1,3-dimethylcyclohexane trans-1,4-dimethylcyclohexane ethylcyclohexane isobutyl isobutyrateb 2,3-dimethylhexane 2-methyl-3-ethylpentane n-octaneb 2,2,3-trimethylpentane 2,2,4-trimethylpentaneb 2,3,3-trimethylpentane 2-ethyl-1-hexanol quinolined cumene o-ethyltoluene 29.00 mesitylene n-propylbenzene 1,2,3-trimethylbenzeneb 1,2,4-trimethylbenzene n-propylcyclohexane 3,3-diethylpentane n-nonane 2,2,3,3-tetramethlpentaneb 1,2,3,4-tetrahydronaphthalene n-butylbenzene p-cymene isobutylbenzeneb n-decane n-nonadecane n-tetradecane n-octadecane sec-butyl chloride methyl isopropyl ether sec-butylamine tert-butylamine 2-methyl-1-butene 2-methyl-2-butene 3-methyl-1-butene methyl isopropyl ketone methyl n-butyrate piperidineb ethyl propyl ether methyl tert-butyl ether 3-methylpyridine 1,5-hexadiene hexanenitriled 3-hexanone ethyl isobutyrate n-propyl propionate 1-hexanol 2-hexanol
exptl pressure linear predicted CNN model predicted exptl temp linear predicted CNN model predicted 29.88 31.58 35.83 45.89 40.55 45.00 49.40 41.00 37.00 33.53 28.70 27.04 26.98 27.77 35.62 34.95 34.65 42.44 35.80 29.00 29.00 29.00 29.00 29.00 30.00 25.76 25.96 26.65 24.57 26.94 25.34 27.83 26.94 45.99 31.67 30.00 29.00 30.86 31.58 34.09 31.90 27.70 26.40 22.60 27.00 35.73 28.49 28.00 30.00 20.82 11.94 15.49 12.53 38.49 38.29 39.48 37.90 33.56 33.56 34.70 38.00 34.28 45.90 33.26 33.85 43.23 33.06 28.82 32.77 30.00 30.20 34.64 33.56
30.45 30.44 30.88 44.58 40.20 43.35 44.37 39.68 41.81 32.91 28.14 26.47 25.80 26.96 35.94 34.85 35.98 40.93 34.46 32.15 30.74 28.89 31.82 30.01 29.93 23.91 23.98 24.28 23.68 27.01 24.69 29.52 27.11 45.11 32.61 32.54 32.54 29.58 31.49 34.85 32.98 27.30 24.96 21.13 28.25 36.28 27.24 30.22 29.33 18.83 15.97 12.76 13.05 38.46 39.38 40.51 36.40 35.62 34.93 35.46 38.71 33.90 43.48 33.79 34.61 43.55 34.35 34.61 33.05 30.36 28.39 36.40 35.12
28.50 32.16 31.98 46.17 39.95 47.39 46.08 39.70 42.11 32.63 28.25 27.34 27.56 28.55 34.44 33.69 34.74 41.55 32.13 30.36 29.73 28.35 30.11 29.46 29.22 26.24 26.30 26.22 24.50 27.39 26.49 28.62 28.16 47.15 31.14 31.50 30.97 29.72 30.05 34.93 31.92 27.36 26.17 21.89 27.65 35.57 27.36 29.20 28.18 19.43 15.46 13.33 12.32 37.26 38.38 40.87 38.45 34.50 33.99 34.28 38.98 34.04 44.28 32.53 35.02 43.39 32.47 32.45 32.74 30.71 29.72 35.25 33.49
530.6 523.1 555.8 695.0 591.8 705.8 697.5 709.1 694.1 569.5 537.3 540.3 530.4 535.3 617.2 617.0 616.3 701.0 687.1 606.1 596.1 591.1 598.0 590.1 609.1 602.0 563.4 567.0 568.8 563.5 544.0 573.5 640.3 782.2 631.1 651.1 640.1 637.4 638.4 664.5 649.1 639.1 610.0 595.6 610.8 720.1 660.5 653.1 650.1 618.5 658.2 692.4 745.3 520.6 464.5 514.3 483.9 465.0 471.0 450.4 553.0 554.5 594.0 500.2 497.1 645.0 507.0 622.1 582.8 553.1 578.0 611.3 586.2
522.5 532.9 553.2 668.8 589.6 701.8 690.7 699.0 700.1 582.5 532.4 532.4 522.5 532.2 620.1 613.3 616.0 715.8 673.9 609.7 605.3 595.8 615.7 598.0 610.0 609.0 557.1 561.2 560.7 559.8 529.3 586.8 629.7 NA 637.3 649.4 644.7 635.1 640.6 654.8 647.1 631.9 613.9 587.1 598.9 710.7 663.5 661.0 651.7 610.7 656.2 689.8 748.8 529.9 467.5 519.7 494.4 472.6 469.0 466.8 554.4 552.8 586.3 493.3 496.7 636.5 510.9 NA 587.1 572.3 576.9 620.6 605.2
522.5 532.9 549.2 687.2 587.2 704.9 700.6 706.5 707.7 580.7 534.0 540.5 527.6 533.1 614.4 608.4 610.6 709.6 683.0 604.1 600.0 592.0 612.6 593.7 606.5 599.5 555.1 558.8 571.1 558.5 536.5 584.4 632.1 NA 633.1 647.2 640.8 629.9 635.4 661.8 646.8 628.9 613.2 596.6 602.0 709.0 660.8 666.5 648.2 618.8 660.7 694.3 742.4 533.8 470.8 521.9 492.8 467.6 463.7 463.7 550.4 556.7 579.9 494.4 500.5 636.9 511.1 NA 579.4 563.6 573.4 611.4 593.9
PREDICTION
OF
CRITICAL TEMPERATURES
AND
PRESSURES
J. Chem. Inf. Comput. Sci., Vol. 38, No. 4, 1998 643
Table 1 (Continued) 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 a
compound name
exptl pressure
linear predicted
CNN model predicted
exptl temp
linear predicted
CNN model predicted
4-methyl-2-pentanol acetonea 2-propanola ethyl vinyl ethera isovaleric acida n-butyl acetatea p-cresola n-toluidinea 1-heptenea methylcyclohexanea 2,2,3-trimethylbutanea o-xylenea 1-octenea m-ethyltoluenea p-diethylbenzenea 1-decenea propionitrilea 4-methylpyridinea n-pentyl formatea
34.25 46.40 53.06 40.17 38.39 30.69 50.83 39.48 27.93 34.26 30.99 36.85 25.17 28.00 27.66 21.40 41.30 46.00 30.84
35.28 50.56 51.03 38.73 37.81 29.55 44.65 41.03 28.73 32.76 29.13 37.41 25.19 30.27 30.01 18.91 47.88 43.11 31.19
32.78 45.65 48.64 38.71 37.95 29.01 50.45 43.43 27.38 31.58 28.69 37.53 24.15 29.22 28.93 19.80 46.39 42.82 28.36
574.4 508.2 508.3 475.1 634.0 579.1 704.6 693.1 537.3 572.2 531.2 630.4 566.6 637.1 658.0 617.0 564.4 646.1 576.0
594.8 471.2 493.7 477.2 635.3 565.8 707.0 704.3 537.4 578.5 527.2 624.0 565.4 642.9 671.6 614.6 560.3 634.0 596.1
584.6 484.5 498.2 476.9 635.6 576.1 707.0 709.3 543.9 578.2 529.8 618.5 573.6 638.1 675.9 620.5 559.6 634.8 582.6
Member of the prediction set. b Member of the cross-validation set. c Outlier of the pressure model. d Outlier of the temperature model.
generation of the model, was used for validation by determining its rms error and comparing it to that of the tset. If the rms error values varied greatly, the model was determined to be invalid. If any of the methods of validation failed, a different model was generated and tested. Neural Networks. Once valid linear models were found using MLR, steps were taken to see if prediction results could be improved by the use of CNN. Typically, superior models can be found using CNNs because they implement nonlinear relationships and because they have more adjustable parameters than the linear models. Previous papers have presented the use of CNNs for QSPR studies in some detail.17 The present work used fully connected, three-layer, feed-forward CNNs which were trained using a quasi-Newton algorithm. The data set was divided into a tset, a cvset, and a pset of 132, 15, and 18 compounds, respectively. These data were then submitted to analysis by networks of various architectures, using the same set of descriptors, to produce the best results with the fewest adjustable parameters. Since the number of input neurons is fixed at the number of descriptors being used and the single output neuron is fixed, investigations of neural network architecture were confined to seeking the best number of hidden layer neurons to use. The cvset was used to locate the optimal set of weights and biases for a given network system and also to avoid overtraining. All computations were performed on a DEC 3000 AXP Model 500 workstation utilizing ADAPT software.3 RESULTS AND DISCUSSION
Critical Temperature Investigations. The first modeling step of this study was the development of a linear model to predict critical temperatures from structural information alone. Two nitrogen-containing compounds (quinoline and hexanenitrile) were determined to be statistical outliers in preliminary modeling, and they were removed from the training set. The eight descriptors comprising the best linear model are shown in Table 2. This model produced an rms error of 9.16 K and a correlation coefficient (R) of 0.993 for the training set compounds. The external prediction set had a rms error of 12.44 K. However, a large percentage of
Table 2. Critical Temperature Model Developed from ADAPT Descriptors Using Multiple Linear Regression explanation
descriptor
coefficients
the dipole moment total charge weighted partial positive surface area simple path 3 connectivity index number of oxygens secondary sp3 carbons cube root gravitation index ∑ ((S.A. * charge of acceptor atoms)/ number of acceptor atoms) average charge on positively charged carbons y intercept
DPOL PPSA 2
13.24 ( 1.39 -0.1025 ( 0.0139
S3P 27.69 ( 2.66 NO -40.68 ( 1.96 2SP3 4.272 ( 0.527 GRAV 3 59.37 ( 2.51 SCAA 2 -7.427 ( 0.192 ACPC CONS
167.2 ( 15.9 73.22
the pset error was associated with one compound, acetone, and when it was removed from the pset the rms error dropped to 9.13 K. Thus, this linear model validates quite well. The linear model described above showed improved results compared to that of a previous study which used the same data set of 165 compounds with an 18 compound pset and a 147 compound tset.2 The previous study of this data set incorporated the experimental boiling point of each compound as a descriptor. That model had an rms error of 8.36 K, R of 0.994, and was validated with an external pset rms error of 8.75 K. In addition to the model using boiling points, another model was constructed using structural information alone. That model had an rms error of 11.51 K and R of 0.989. However, the overall rms error for the pset was not reported. Neither of the original models contained any statistical outliers. The eight descriptors in the present linear model span the range of descriptor types. Three are topological descriptors (S3P, NO, 2SP3), two are electronic descriptors (DPOL, ACPC), one is a geometric descriptor (GRAV 3), one is a charged partial surface area descriptor (PPSA 2), and one is a hydrogen-bonding descriptor (SCAA 2). The three most important descriptors of the set, by their F-values, are S3P, PPSA 2, and SCAA 2. These eight descriptors, taken as a group, successfully encode the compounds’ structures. Comparison of the two structural models proves to be valuable. The models contain two descriptors in common:
644 J. Chem. Inf. Comput. Sci., Vol. 38, No. 4, 1998
TURNER
ET AL.
Table 3. Critical Pressure Model Developed from ADAPT Descriptors Using Multiple Linear Regression
Figure 1. Critical temperature calculated vs observed plot for the training, cross-validation, and prediction sets using an eightdescriptor computational neural network model.
the number of oxygens and the average charge on positively charged carbons. This indicates the possible importance of the relationship between these descriptors and critical temperature. In addition, both models contained descriptors of some type of path connectivity and charge normalized by surface area. Regardless of the similarities, the results of the two models greatly differ with a 36% decrease in training set rms error between the new and the old model. It is believed that the best explanation for the decreased error comes from the presence of new descriptors. Three descriptors not previously available are used here: 2SP3, GRAV 3, and SCAA 2. These new descriptors evidently capture structural information in a new way that was not possible with the descriptors available to the previous study. Having shown that a linear model could be constructed from structural descriptors alone, the next step involved using a CNN to develop a nonlinear model based on the same set of descriptors. The eight descriptors were tested with several CNN architectures, and it was determined that a 8-4-1 architecture (8 input neurons for the eight descriptors, 4 hidden neurons, and 1 output neuron for a total of 41 adjustable parameters) yielded the best results, considering the accuracy of prediction and the number of adjustable parameters. The CNN calculations were successful in improving the results. The observed rms errors were 7.33 K for the tset, 7.76 K for the cvset, and 9.85 K for the pset. However, once again, acetone contributed a large percentage of the error in the prediction set. Once acetone was removed the rms error of the pset was reduced to 8.34 K. The 7.33 K rms error for the training set corresponds to a 1.3% relative error at the mean of the critical temperature values (573.5 K). The calculated versus experimental critical temperatures of the CNN model are shown in Figure 1. The plot shows that the model is of high quality throughout the range of the critical temperature values. Examination of the rms errors for both the linear and CNN models illustrates the improved results from the previous study. The present study, using structural information alone,
explanation
descriptor
the charge on the most negative atom partial positive surface area weighted partial negative surface area the Wiener number number of paths length 3 second major moment (width) first major moment of inertia/ second major moment of inertia sum of surface area of acceptor atoms/number of acceptor atoms y intercept
QNEG 1
coefficients -14.57 ( 1.83
PPSA 1 WNSA 2
-0.1244 ( 0.0059 0.4344 ( 0.052
ALLP 5 N3P 5 GEOM 2 MOMI 4
0.03641 ( 0.0047 0.7220 ( 0.087 -8.536 ( 1.03 10.76 ( 1.26
SAAA 2
0.1201 ( 0.021
CONS
53.46
was able to achieve results better than those determined previously using a physical property as a descriptor. In addition, another indication of excellent results is that the errors are very near the predicted error of 8 K determined by the reported experimental values. Critical Pressure Investigations. The next step was to build a model to predict the critical pressure. This together with the critical temperature would allow the critical point to be determined. The same data set of compounds was employed with the same training, cross-validation, and prediction sets. As before, an eight-descriptor model proved to be the most effective for prediction. However, phenol proved to be an outlier in this phase of the study. The eight descriptors comprising the final model are shown in Table 3. This linear model produced a rms error of 2.07 atm and R of 0.964. The pset rms error was 2.80 atm. The majority of the model’s error was found in the prediction of compounds with high critical pressure values. The eight descriptors chosen for modeling of critical pressure span the range of descriptor types available for use. Two of the descriptors are topological (ALLP 5, N3P 5), 1 is electronic (QNEG 1), 2 are geometric (GEOM 2, MOMI 4), 2 are charged partial surface area descriptors (PPSA 1, WNSA 2), and 1 is a hydrogen-bonding descriptor (SAAA 2). Only one descriptor is found in both the critical temperature and critical pressure models, the hydrogenbonding descriptor SAAA 2. The three most important descriptors of the set, by their F-values, are PPSA 1, QNEG 1, and ALLP 5. This set of eight descriptors, taken as a group, successfully encodes the structures of the compounds. Once the linear model was developed, a CNN was used to improve the overall results. The best CNN architecture was determined to be 8-5-1 (51 adjustable parameters). The best CNN model found gave a tset rms error of 1.51 atm, a cvset rms error of 1.35, and an external pset error of 2.39. The 1.51 atm rms error for the training set corresponds to a 4.4% relative error at the mean of the critical pressure values (33.7 atm). Although the pset error still remains well above that of the tset, it can be seen in Figure 2 that the general predictive capability of the CNN model is quite good. The predictions in the higher pressure ranges remain good. This plot shows more scatter than that for the critical temperatures; for this set of compounds, the critical pressures are more difficult to model.
PREDICTION
OF
CRITICAL TEMPERATURES
AND
PRESSURES
J. Chem. Inf. Comput. Sci., Vol. 38, No. 4, 1998 645 REFERENCES AND NOTES
Figure 2. Critical pressure calculated vs observed plot for the training, cross-validation, and prediction sets using an eightdescriptor computational neural network model.
CONCLUSIONS
This study shows that it is possible to develop good quantitative structure-property relationships for the prediction of critical temperature and pressure of structurally diverse organic compounds using structural information alone. The models provide some insight into what structural features are related to these specific physical properties. Linear models chosen using the genetic algorithm and simulated annealing feature selection methods were able to predict these two critical properties. Additionally, using CNNs to build nonlinear models based on these same sets of descriptors produced even better models with good predictive ability. ACKNOWLEDGMENT
We thank Matthew Wessel, Jon Sutter, Brooke Mitchell, Stephen Johnson, and Heidi Engelhardt for continued assistance throughout this project. This work was supported in part by the Design Institute for Physical Property Data (DIPPR) Project 931: Data Prediction Methods. Support was also received from the National Science Foundation REU Program.
(1) Bea, H.-K.; Lee, S.-Y. A Method for the Prediction of Critical Temperature and Pressure of Pure Fluids. Fluid Phase Equilibria 1991, 66, 225-232. (2) Egolf, L. M.; Wessel, M. D.; Jurs, P. C. Prediction of Boiling Points and Critical Temperatures of Industrially Important Organic Compounds from Molecular Structure. J. Chem. Inf. Comput. Sci. 1994, 34, 947-956. (3) Stuper, A. J.; Brugger, W. E.; Jurs, P. C. Computer-Assisted Studies of Chemical Structure and Biological Function; Wiley-Interscience: New York, 1979; p 83. (4) Johnson, S. R.; Jurs, P. C. Acute Mammalian Structure for a Diverse Set of Substituted Anilines Using Regression Analysis and Computational Neural Networks. In Computer-Assisted Lead Finding and Optimization; van de Waterbeemd, H., Testa, B., Folkers, G., Eds.; Verlag Helvetica Chimica Acta: Basel, 1997. (5) Mitchell, B. E.; Jurs, P. C. Prediction of Autoignition Temperatures of Organic Compounds from Molecular Structure. J. Chem. Inf. Comput. Sci. 1997, 37, 538-547. (6) Engelhardt, H. L.; Jurs, P. C. Prediction of Supercritical Carbon Dioxide Solubility of Organic Compounds from Molecular Structure. J. Chem. Inf. Comput. Sci. 1997, 37, 478-484. (7) Wessel, M. D. Computer-Assisted DeVelopment of QuantitatiVe Structure-Property Relationships and Design of Feature Selection Routines. Ph.D. Dissertation, Pennsylvania State University, University Park, PA, 1996. (8) Sutter, J. M.; Jurs, P. C. Selection of Molecular Structure Descriptors for Quantitative Structure-Activity Relationships. In Adaption of Simulated Annealing to Chemical Problems; Kalivas, J. H., Ed.; Elsevier Science Publishers B. V.: Amsterdam, 1995. (9) Design Institute for Physical Property Data (DIPPR). Physical and Thermodynamic Properties of Pure Chemicals: Data Compilation; Daubert, T. E., Danner R. P., Eds.; Hemisphere Publishing: New York, 1989: Vol. 1-4. (10) Hypercube, Inc., Waterloo, ON. (11) Stewart, J. P. P. MOPAC 6.0, Quantum Chemistry Program Exchange; Indiana University, Bloomington, IN, Program 455. (12) Stewart, J. P. P. MOPAC, A Seimemperical Molecular Orbital Program. J. Comput.-Aided Mol. Des. 1990, 4, 1. (13) Stanton, D. T.; Jurs, P. C. Development and Use of Charged Partial Surface Area Structural Descriptors for Quantitative StructureProperty Relationship Studies. Anal. Chem. 1990, 62, 2323. (14) Russell, C. J.; Dixon, S. L.; Jurs, P. C. Computer Assisted Study of the Relationship Between Molecular Structure and Henry’s Law Constant. Anal. Chem. 1992, 64, 1350. (15) Belsey, D. A.; Kuh, E.; Welsch R. E. Regression Diagnostics; Wiley: New York, 1980. (16) Stanton D. T.; Jurs, P. C. Computer Assisted Prediction of Normal Boiling Points of Furans, Tetrahydrofurans, and Thiophenes. J. Chem. Inf. Comput. Sci. 1991, 31, 301. (17) Xu, L.; Ball, J. W.; Dixon, S. L.; Jurs, P. C. Quantitative Structureroperty Relationships for Toxicity of Phenols Using Regression Analysis and Computational Neural Networks. EnViron. Toxical. Chem. 1994, 13, 841.
CI9800054
