A Modified 3D-QSAR Model Based on Ideal Point Method and Its Application in the Molecular Modification of Plasticizers with Flame Retardancy and Eco-Friendliness

The addition of plasticizers makes plastics flammable, and thus, poses a potential risk to the environment. In previous researches, plasticizers with flame retardancy had been synthesized, but their eco-friendliness had not been tested or described. Thus, in this paper, eco-friendliness plasticizers with flame retardancy were designed based on phthalic acid esters (PAEs), which are known as common plasticizers and major plastic additives. For a comprehensive analysis, such as flammability, biotoxicity, and enrichment effects, 17 PAEs’ comprehensive evaluation values were calculated based on the ideal point method. Further, a multi-effect three-dimensional quantitative structure-activity relationship (3D-QSAR) model of PAEs’ flammability, biotoxicity and enrichment effects was constructed. Thus, 18 dimethyl phthalate (DMP) derivatives and 20 diallyl phthalate (DAP) derivatives were designed based on three-dimensional contour maps. Through evaluation of eco-friendliness and flammability, six eco-friendly PAE derivatives with flame retardancy were screened out. Based on contour maps analysis, it was confirmed that the introduction of large groups and hydrophobic groups was beneficial to the simultaneous improvement of PAEs’ comprehensive effects, and multiple effects. In addition, the group properties were correlated significantly with improved degrees of the comprehensive effects of corresponding PAE derivatives, confirming the feasibility of the comprehensive evaluation method and modified scheme.


Introduction
Numerous properties of plastic materials, lightweight features, corrosion resistance ability, bright color, transparency, ease to process, etc., have made them greatly applicable in many fields. Furthermore, the versatility, low cost, and abundance of plastic-based products also influence their daily applicability [1]. Due to being constructed by organic molecules, plastics are rich with ignition properties, which contribute to a potential fire-risk [2]. Due to the highly flammable property of plastic materials, they have a restricted use in high-temperature areas. Therefore, to avoid the significant loss of lives and properties caused by fire disaster, it is necessary to reduce the flammability of plastics through appropriate flame retardant treatment (such as covering the surface with a flame retardant coating and adding flame retardants) [1,3]. However, the addition of flame retardants is one of the prime methods to eliminate the flammability of plastics. So far, many flame retardants have been used to improve the flame retardant performance of plastics, which includes inorganic, halogenated, In summary, plastics and its main additive PAEs are consistent in exhibiting several adverse effects in the environment and eco-system that include flammability, biotoxicity, and enrichment. On the other hand, the application of many flame retardant additives and their byproducts cause serious health hazards and environmental issues, which limit their industrial application [14,24]. Due to the same reason, development of low toxic and eco-friendly flame retardant technologies has become a very important research direction. In this paper, PAEs were used as the prime research objects to obtain information about the flammability property, biotoxicity, and enrichment effects of plastics, according to which new eco-friendly PAE-based plasticizers with flame retardancy were developed. In order to evaluate the molecules comprehensively, it is necessary to treat a variety of indicators with different dimensions using mathematical methods or equations. For example, Vahabi et al. [25,26] defined a simple universal dimensionless index for the first time, known as the flame retardancy index (FRI), which was useful for the comparative evaluation of the flame retardancy performance of thermoplastic systems regardless of the types of polymers and additives used. Li et al. [27] evaluated molecules' multi-effect comprehensively via a simplified formula method and standard deviation score method. In this paper, dimensionless processed activity data of the flammability, biological (fish) toxicity and enrichment effects of 17 PAEs were evaluated via mathematical methods (ideal point method), which further combined to calculate the multi-effect comprehensive evaluation values of 17 PAEs. Then, a three-dimensional quantitative structure-activity relationship (3D-QSAR) model of PAEs' flammability, biotoxicity, and enrichment multi-effect was constructed by taking the PAEs' flammability, biotoxicity and enrichment multi-effect comprehensive evaluation values as the dependent variables and the molecular structure parameters as the independent variables. To obtain PAE derivatives with low biotoxicity, low enrichment, and flammability, the target molecules dimethyl phthalate (DMP) and diallyl phthalate (DAP) were modified by single and double substituent groups through the three-dimensional contour map of the multi-effect model. In addition, through the prediction and evaluation of PAEs' multi-effect 3D-QSAR model, flammability, biotoxicity, enrichment, degradability, and mobility single-effect 3D-QSASR models, it was confirmed that the designed PAE derivatives were eco-friendly molecules with the better flame retardant property. At the same time, the stability and insulation of PAE derivatives were also evaluated.

Sources of Data
Most PAEs are used as plasticizers in the processing of polyvinyl chloride (PVC) [27]. Referring to Li et al. [27], 17 kinds of PAEs from PVC plastics were selected as representatives, which consist of higher molecular weight molecules and lower molecular weight molecules, long side-chain molecules and short side-chain molecules, different structural molecules including ring, linear chain or branched chain. The flammability property of 17 kinds of PAEs was characterized by LOI. Greater LOI value indicates the requirement of a higher amount of ultimate oxygen concentration for the combustion of material, and subsequently confirms the better flame retardant capability [18]. The flame retardant property of polymer materials is related to the molecular structure of the polymer itself. Empirical Formulas (1) and (2), which are applicable for halogen-free materials to calculate the ratio between the residue amount and flame retardancy during thermal decomposition, are commonly used to calculate the LOI value [28]. Following the same rule, here in this research work, the LOI data of PAEs were obtained from the calculations of Formulas (1) and (2) ( Table 1). The biotoxicity data of 17 PAEs was obtained through the Estimation Program Interface (EPIWEB) 4.1 (Copyright© 2020-2012, United States Environmental Protection Agency for EPI Suite TM and all component programs except BioHCWIN and KOAWIN, U.S.) software and employed as the 50% lethal concentration of PAEs to fish (lethal concentration of 50%, LC 50 ) ( Table 1). Bio-concentration factors (logBCF) of 17 PAEs were predicted by EPIWEB software to characterize the molecular enrichment (Table 1). Half-life (HL) data of 17 PAEs in rivers were introduced from the EPIWEB software; the logarithm value of the same (logHL) was used as an indicator of degradation abilities of the PAEs (Table 1). EPIWEB software was used to obtain PAEs' experimental or predictive values for octanol/air partition coefficient (logK OA ) to represent their migration (Table 1). In addition, the hydrophobicity parameter octanol/water partition coefficient (logP) of PAEs derivatives was calculated by software Chembiodraw 12.0.
where CR refers to the amount of residue when the substance is heated at 850 • C.
In the formula, (CFT) i is the contribution index of the ith functional group to the residue amount of a substance. The CFT value of the functional group Polymers 2020, 12, x FOR PEER REVIEW 4 of 21 characterize the molecular enrichment (Table 1). Half-life (HL) data of 17 PAEs in rivers were introduced from the EPIWEB software; the logarithm value of the same (logHL) was used as an indicator of degradation abilities of the PAEs (Table 1). EPIWEB software was used to obtain PAEs' experimental or predictive values for octanol/air partition coefficient (logKOA) to represent their migration (Table 1). In addition, the hydrophobicity parameter octanol/water partition coefficient (logP) of PAEs derivatives was calculated by software Chembiodraw 12.0.
where CR refers to the amount of residue when the substance is heated at 850 °C.
In the formula, (CFT)i is the contribution index of the i th functional group to the residue amount of a substance. The CFT value of the functional group contained in PAEs is 2, while the CFT values of other halogen-free functional groups are 0. M refers to the molecular weight of monomers in the polymer, and here refers to the relative molecular weight of PAEs. To eliminate the influence of the index dimension on molecular evaluation, first, the performance indicators of PAEs' flammability, biotoxicity, and enrichment should be processed dimensionless using Formula (3), (4) via normalization method, which is the superior way for the comprehensive evaluation of PAEs. The normalized index data Cij, is processed by the ideal point method (Formula (5)) [29]; further, the comprehensive evaluation value Z of PAEs' flammability, biotoxicity and enrichment was calculated by weighting them together, with the weight ratio of PAEs' flammability, biotoxicity, and enrichment at 40%: 30%: 30%. The above step should be performed before constructing a 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment multi-effect, so as to evaluate the PAEs comprehensively.
contained in PAEs is 2, while the CFT values of other halogen-free functional groups are 0. M refers to the molecular weight of monomers in the polymer, and here refers to the relative molecular weight of PAEs. To eliminate the influence of the index dimension on molecular evaluation, first, the performance indicators of PAEs' flammability, biotoxicity, and enrichment should be processed dimensionless using Formulas (3) and (4) via normalization method, which is the superior way for the comprehensive evaluation of PAEs. The normalized index data C ij , is processed by the ideal point method (Formula (5)) [29]; further, the comprehensive evaluation value Z of PAEs' flammability, biotoxicity and enrichment was calculated by weighting them together, with the weight ratio of PAEs' flammability, biotoxicity, and enrichment at 40%:30%:30%. The above step should be performed before constructing a 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment multi-effect, so as to evaluate the PAEs comprehensively.
Polymers 2020, 12, 1942 5 of 18 where C ij refers to the ith dimensionless parameter of the jth index; maxX j and minX j are the maximum and minimum values of the jth index, respectively; for the positive indexes, it is applicable to Formula (3); for the inverse indexes, Formula (4) is applicable.
where λ j is the weight of jth index, among which λ 1 = 40%, λ 2 = 30%, λ 3 = 30%; C j * is the ideal point value of the jth index. After normalization of PAEs' effect index parameters, 0 ≤ C ij ≤ 1, and all of them were positive indexes; therefore, the value of C j * in the equation is 1. Z i is called Euclidian distance, namely the comprehensive evaluation value of the ith PAEs' flammability, biotoxicity, and enrichment.

Construction of PAEs' Flammability, Biotoxicity, and Enrichment Multi-Effect 3D-QSAR Model
The three-dimensional structures of PAEs were drawn in Sybyl-X 2.0; molecular geometry was mechanically optimized to obtain the low-energy conformation of each molecule using the Powell method followed by the Conjugate Gradient method while the Minimize module was employed. The molecular charge was set as Gasteiger-Huckle. Molecular force field Tripos was exerted on molecules, and the energy convergence standard was set as 0.005 kcal·mol −1 with 10,000 iterations.
The most active molecule was selected as the template molecule, while the Align Database module was used to superposition the skeleton of the compounds. During the calculation of the parameters of the Comparative Molecular Similarity Index Analysis (CoMSIA) field, the types of molecular force fields were selected as steric field (S), hydrophobic field (H), electrostatic field (E), hydrogen bond donor field (D), and hydrogen bond acceptor field (A), to gain the molecular field parameters. The dielectric constant was related to the distance when the threshold value was set as 125.4 kJ/mol. The remaining parameters adopted the default values of the system. CoMSIA field parameters could be calculated automatically via Sybyl-X 2.0 by entering PAEs' Eigenvalues of the training set. When the Partial Least Square (PLS) method was used to analyze the relationship between the molecular structure and biological activity [30], the Leave-One-Out (LOO) method was used to cross-validate the training set compounds initially [31], which is required to estimate the cross-validation coefficient q 2 and the optimal number of principal components n. Thereafter, "No Validation" was used for regression analysis while r 2 , standard error of estimate (SEE), and Fischer's value (F) were calculated to complete the construction of the CoMSIA model [32]. In addition, the CoMSIA model needs to be verified externally by external test sets. The interactive validation coefficient r 2 pred and standard deviation standard error of predict (SEP) of the test sets were used for confirming the external prediction ability of the model [33][34][35].

Evaluation of Eco-Friendliness, Stability, and Insulation of PAE Derivatives
In this paper, the single-effect 3D-QSAR models of PAEs' biotoxicity, enrichment, degradability, and mobility were established, respectively, for testing the eco-friendliness of PAE derivatives. The single-effect CoMSIA models of PAEs' biotoxicity, enrichment, degradability, and mobility were constructed by considering PAEs' structure parameters as independent variables and PAEs' logLC 50 , logBCF, logHL and logK OA values as dependent variables. To screen the eco-friendliness of PAE derivatives, logLC 50 , logBCF, logHL and logK OA values of the designed PAE derivatives were predicted using the above model.
Using software Gaussian 09, the structure of PAEs and its derivatives was optimized at the level of B3LYP/6-31G* base group via density functional theory (DFT). The positive frequency value (representing stability) and energy gap value (representing insulation) were calculated subsequently, while the results of the same were viewed by GaussView 5.0 and Multiwfn. When the positive frequency (Freq.) of PAE derivative exhibits a value greater than zero, it indicates the molecular structure of the designed PAE derivative is stable and can exist in a stable form in the environment [36]. The energy difference between the highest occupied molecular orbital (HOMO) and the lowest unoccupied molecular orbital (LUMO) is called Energy gap, which refers to the transition ability of the electrons from an occupied orbital to vacant orbital; this is an important parameter reflecting the electrical conductivity and luminescence properties of materials. A greater energy gap value indicates a weaker conductivity of the molecule [37].

Construction and Verification of Multi-Effect 3D-QSAR Model of PAEs' Flammability, Biotoxicity, and Enrichment
In this paper, a 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment multi-effect was constructed by considering the molecular structure parameters of PAEs as independent variables and the comprehensive evaluation value Z of PAEs' flammability, biotoxicity and enrichment multi-effect as dependent variables.
The training set is used to construct the PAEs' multi-effect 3D-QSAR model, and the test set is used to verify its accuracy. According to the ratio of 3:1, 13 PAEs were randomly selected as the training sets, and the remaining four PAEs were used as the test sets. DMEP was used as the template molecule to build the multi-effect 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment; the results are shown in Table 2. The evaluation parameters and molecular force field contribution rate of PAEs' multi-effect CoMSIA model are shown in Table 3.  As per the results in Table 3, the optimal principal component n of PAEs' multi-effect CoMSIA model was 10, while the cross-validation coefficient q 2 was 0.616 (generally considered that when the value of q 2 > 0.5, the model has the credible predictive ability [38,39]). The non-cross-validation coefficient r 2 was 1.000 (>0.9), SEE was 0.008, and F was 507.314, indicating the greater fitting ability and robustness of the model [40,41]. The external test coefficient r 2 pred was 0.66 (>0.6), and SEP of the test set was 0.16, which specified the model also had well predictive ability [42,43] and the constructed CoMSIA model met the requirements [44]. The model was used to predict the comprehensive evaluation value Z of PAE derivatives' flammability, biotoxicity, and enrichment multi-effect. In addition, from the contribution rates of the molecular force fields in the CoMSIA model in Table 3, it could be seen that the contribution rates of S, E, H, D and A were 22.9%, 14.0%, 44.9%, 0.00%, and 11.2% separately, which indicated the spatial effect, electrical distribution, and hydrophobicity of the groups had a large impact on PAEs' flammability, biotoxicity, and enrichment effects. In contrast, the effect of the hydrogen bond donor field was negligible.
In order for further verification, the flame retardancy of PAE derivatives, and a 3D-QSAR model of PAEs' flammability single-effect was also established in this paper. In addition, single-effect 3D-QSAR models of PAEs' biotoxicity, enrichment, degradability, and mobility were established to evaluate the eco-friendliness of PAE derivatives. In accordance with the above steps, the 3D-QSAR model of PAEs' flammability single-effect was constructed and tested by selecting 13 PAEs arbitrarily as the training set, while the remaining four PAEs were considered as the test set; DMP used as the template molecule. Finally, the CoMSIA model that met the requirements was obtained. The optimal principal component n of the CoMSIA model was 10, the cross-validation coefficient q 2 was 0.580 (>0.5), the non-cross-validation coefficient r 2 was 1.000 (r 2 > 0.9), SEE was 0.010, and F value was 9275.112, which indicated the model had good fitting ability and robustness. The external test coefficient r 2 pred was 0.901 (>0.6), and SEP was 0.346, indicating that the model had a good predictive ability. So, the constructed CoMSIA model could be used to predict the flammability of PAEs derivatives.
In addition, PAEs' single-effect CoMSIA models of biotoxicity, enrichment, degradability, and mobility were also established in the research sequence, which were used to test the eco-friendliness of the designed molecules. The 3D-QSAR model of PAEs' biotoxicity single-effect was established by considering PAEs' molecular structure parameters as independent variables and logLC 50 of PAEs on fish as the dependent variables. The 3D-QSAR model of PAEs' enrichment single-effect was established by considering PAEs' molecular structures as independent variables and logBCF of PAEs as the dependent variables. The 3D-QSAR model of PAEs' degradability single-effect was established, by considering PAEs' molecular structure information as the independent variable and the degradable activity data (logHL) of PAEs as the dependent variable. The 3D-QSAR model of PAEs' mobility single-effect was established by considering PAEs' molecular structure parameters as the independent variables and logK OA of PAEs as the dependent variables. By evaluating the above parameters, it may be granted that the constructed models met the requirements and had good predictive abilities ( Table 4). The contribution rates of the models' molecular force fields are shown in Table 5. The experimental values, predicted values, and relative errors of the molecules in the model training set and test set are given in Table 6.  Table 5. Molecular force field contribution rates of PAEs' flammability, biotoxicity, enrichment, degradability, and mobility models.

Molecular Modification and Evaluation Based on PAEs' 3D-QSAR Model
To obtain flame retardant PAE derivatives with low biotoxicity and low enrichment, the CoMSIA model contour maps were generated by a multi-effect 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment; molecular modification schemes for single and double substitutions of target molecules were formulated, by selecting DMP and DAP as target molecules. The three-dimensional contour maps explain the interaction of small molecules, the relationship between the molecules, and the receptors visibilities, which provide a basis for modifying compounds; the colorful blocks of the contour maps are beneficial to enhance molecular activity by default. Taking DMP as an example, the three-dimensional contour maps of the CoMSIA model were analyzed (Figure 1). The yellow block of the steric field indicated that the introduction of smaller groups in this area can improve the activity of DMP derivatives. Therefore, the introduction of bulky groups here was conducive to reduce the comprehensive effect value of the molecule, since PAEs' comprehensive evaluation value was an inverse index. In the electrostatic field, the introduction of positively charged groups in the blue region, and the introduction of negatively charged groups in the red region are beneficial to improve the activity of the molecule. In the hydrophobic field, the introduction of hydrophilic groups in the white area helps to improve the activity of the molecule. In order to reduce the comprehensive evaluation value of DMP's flammability, biotoxicity, and enrichment effects, hydrophobic groups were introduced into the molecule. In the hydrogen bond receptor field, the purple-colored block indicates that increasing the number of hydrogen bond receptors is beneficial to improve the activity of DMP, while the red areas indicate the increased number of hydrogen bond donors reduces the activity of molecules.

Molecular Modification and Evaluation Based on PAEs' 3D-QSAR Model
To obtain flame retardant PAE derivatives with low biotoxicity and low enrichment, the CoMSIA model contour maps were generated by a multi-effect 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment; molecular modification schemes for single and double substitutions of target molecules were formulated, by selecting DMP and DAP as target molecules. The three-dimensional contour maps explain the interaction of small molecules, the relationship between the molecules, and the receptors visibilities, which provide a basis for modifying compounds; the colorful blocks of the contour maps are beneficial to enhance molecular activity by default. Taking DMP as an example, the three-dimensional contour maps of the CoMSIA model were analyzed (Figure 1). The yellow block of the steric field indicated that the introduction of smaller groups in this area can improve the activity of DMP derivatives. Therefore, the introduction of bulky groups here was conducive to reduce the comprehensive effect value of the molecule, since PAEs' comprehensive evaluation value was an inverse index. In the electrostatic field, the introduction of positively charged groups in the blue region, and the introduction of negatively charged groups in the red region are beneficial to improve the activity of the molecule. In the hydrophobic field, the introduction of hydrophilic groups in the white area helps to improve the activity of the molecule. In order to reduce the comprehensive evaluation value of DMP's flammability, biotoxicity, and enrichment effects, hydrophobic groups were introduced into the molecule. In the hydrogen bond receptor field, the purple-colored block indicates that increasing the number of hydrogen bond receptors is beneficial to improve the activity of DMP, while the red areas indicate the increased number of hydrogen bond donors reduces the activity of molecules. During the modification of the target molecules, only the steric field and hydrophobic field with the highest contribution rates in the molecular field were considered. The three-dimensional contour maps showed that the substitution sites are surrounded by the steric field and the hydrophobic field, which provided a reference for molecular modification of PAEs. It could be seen in Figure 1A that the yellow regions were distributed near the H atom connected to the 2-C atom and the 2′-C atom in the DMP and the white region in Figure 1C covered the side chain on the benzene ring of DMP, which indicated the introduction of bulky groups or hydrophobic groups at these sites could be capable of reducing the comprehensive evaluation value of PAEs' flammability, biotoxicity, and enrichment multi-effect. When the molecules were modified, it was noteworthy that the basic skeleton part of PAEs should be retained, as well as retain their performance as a plasticizer. Based on the contour map information of the steric field and the hydrophobic field, the substitution sites of the target molecules DMP and DAP are shown in Figure 2. During the modification of the target molecules, only the steric field and hydrophobic field with the highest contribution rates in the molecular field were considered. The three-dimensional contour maps showed that the substitution sites are surrounded by the steric field and the hydrophobic field, which provided a reference for molecular modification of PAEs. It could be seen in Figure 1A that the yellow regions were distributed near the H atom connected to the 2-C atom and the 2 -C atom in the DMP and the white region in Figure 1C covered the side chain on the benzene ring of DMP, which indicated the introduction of bulky groups or hydrophobic groups at these sites could be capable of reducing the comprehensive evaluation value of PAEs' flammability, biotoxicity, and enrichment multi-effect. When the molecules were modified, it was noteworthy that the basic skeleton part of PAEs should be retained, as well as retain their performance as a plasticizer. Based on the contour map information of the steric field and the hydrophobic field, the substitution sites of the target molecules DMP and DAP are shown in Figure 2. The comprehensive evaluation value Z of 38 PAE derivatives was predicted using the 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment multi-effect. LOI, logLC50, logBCF of PAE derivatives were predicted by PAEs' flammability, biotoxicity, and enrichment single-effect 3D-QSAR models, respectively. Finally, by comparing with the target molecules, it was found that only 22 PAE derivatives had a lower comprehensive evaluation value, improved flammability, decreased biotoxicity, and decreased enrichment. The specific prediction results were shown in Table 7.
It could be seen from Table 7 that the 22 kinds of PAE derivatives had an increasing rate of LOI value of 0.74-9.12%, expressing their flame retardant properties were partially improved. The logLC50 values of PAE derivatives were reduced by one to three orders of magnitude in comparison to the target molecules, with a reduced rate of 9.50-338.22%. As for enrichment effect, the reduced degrees of logBCF of PAE derivatives were 10.79-90.99%. The logBCF values were all less than 3.30, indicating that the molecules would not be accumulated in the organisms [45]. When the LOI value remains 22-27%, molecules represent good flame retardant properties, and they are not easy to burn [28]. In order to obtain PAE derivatives with an obvious modification on the decreased flammability effect, seven PAE derivatives expressing flammability improvement rate >5% were selected for the subsequent analysis, such as DAP-2-CH2NO2, DAP-1-NO2-2-CH2C6H5, DAP-1-NO2-2-CH2CH3, DAP-1-NO2-2-CH2NO2, DAP-1-NO2-2-NO2, DAP-1-NO2-2-OCH3, and DAP-2-CH=CH2-1-NO2. Eight kinds of hydrophobic and bulky groups, such as -CH 3 , -CH 2 CH 3 , -CH 2 C 6 H 5 , -NO 2 , -CH 2 NO 2 , -SH, -OCH 3 , -CH=CH 2 , were selected to perform single and double position modifications of target molecules. Further, a total of 18 DMP derivatives and 20 DAP derivatives were designed.
The comprehensive evaluation value Z of 38 PAE derivatives was predicted using the 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment multi-effect. LOI, logLC 50 , logBCF of PAE derivatives were predicted by PAEs' flammability, biotoxicity, and enrichment single-effect 3D-QSAR models, respectively. Finally, by comparing with the target molecules, it was found that only 22 PAE derivatives had a lower comprehensive evaluation value, improved flammability, decreased biotoxicity, and decreased enrichment. The specific prediction results were shown in Table 7.
It could be seen from Table 7 that the 22 kinds of PAE derivatives had an increasing rate of LOI value of 0.74-9.12%, expressing their flame retardant properties were partially improved. The logLC 50 values of PAE derivatives were reduced by one to three orders of magnitude in comparison to the target molecules, with a reduced rate of 9.50-338.22%. As for enrichment effect, the reduced degrees of logBCF of PAE derivatives were 10.79-90.99%. The logBCF values were all less than 3.30, indicating that the molecules would not be accumulated in the organisms [45]. When the LOI value remains 22-27%, molecules represent good flame retardant properties, and they are not easy to burn [28]. In order to obtain PAE derivatives with an obvious modification on the decreased flammability effect, seven PAE derivatives expressing flammability improvement rate >5% were selected for the subsequent analysis, such as DAP-2-CH 2 NO 2 , DAP-1-NO 2 -2-CH 2 C 6 H 5 , DAP-1-NO 2 -2-CH 2 CH 3 , DAP-1-NO 2 -2-CH 2 NO 2 , DAP-1-NO 2 -2-NO 2 , DAP-1-NO 2 -2-OCH 3 , and DAP-2-CH=CH 2 -1-NO 2 .

Evaluation of Eco-Friendliness, Stability, and Insulation of PAE Derivatives
The eco-friendliness (i.e., degradability, long-distance migration) and other properties (i.e., stability and insulation) of the PAE derivatives screened in Section 3.2 were further evaluated. We calculated logHL (representing degradability), logK OA (representing long-distance mobility), Freq. (representing stability) and energy gap (representing insulation) of PAEs and their derivatives, and the results are given in Table 8. The predictive results of PAEs' degradability 3D-QSAR model showed that logHL values of all PAE derivatives were decreased, which explains the enhancement in their degradability. The logHL value of DAP-2-CH=CH 2 -1-NO 2 was 11.84% lower than that of DAP. According to the predicted results of PAEs' mobility 3D-QSAR model, although the logK OA values of PAEs derivatives were reduced (except DAP-2-CH 2 NO 2 ), all the logK OA values of derivatives situated in a range of 6.5-10 with the same molecular mobility level, while all of them were semi-volatile substances [46]. In addition, all the positive frequencies Freq. of PAE derivatives were greater than 0, indicating that the designed molecules could exist stably in the environment. The insulativity of PAEs derivatives was reduced to some extent but within an acceptable range. It was confirmed that the six PAEs derivatives, DAP-1-NO 2 -2-CH 2 C 6 H 5 , DAP-1-NO 2 -2-CH 2 CH 3 , DAP-1-NO 2 -2-CH 2 NO 2 , DAP-1-NO 2 -2-NO 2 , DAP-1-NO 2 -2-OCH 3 , and DAP-2-CH=CH 2 -1-NO 2 , were eco-friendly plasticizers with flame retardancy.

Mechanism Analysis of PAE Derivatives' Improved Multiple Effects Based on the Three-Dimension Contour Maps
In this paper, DMP was taken as an example to compare three-dimension contour maps of PAEs' multi-effect model with flammability, biotoxicity, and enrichment single-effect models, respectively, so as to qualitatively analyze the mechanism of PAE derivatives' improved flammability, biotoxicity, and enrichment comprehensive effect.
As can be seen from Table 5, the steric field and hydrophobic field had the maximum contribution to the force field contribution rate of PAEs' multi-effect CoMSIA model and PAEs' flammability, biotoxicity, enrichment single-effect models. The comprehensive evaluation value Z and the activity data of the multi-effect CoMSIA model were as small as possible, while the modified PAE derivatives showed reduced comprehensive evaluation value Z by changing the structure of the molecules. The characteristic values LOI of flammability and logLC 50 of biotoxicity, both were positive indicators, while the larger values indicated the better results. The smaller values of logBCF indicated lower enrichment of molecules in organisms and thus represented the smaller environmental effects.
In Figure 3, for the steric field, a three-dimension contour map in the multi-effect model was shown in yellow, indicating that the introduction of small groups in the molecule would increase the comprehensive evaluation value. Therefore, in order to reduce the comprehensive evaluation value, larger groups should be introduced into the overlapping part of DMP with three-dimension contour maps during modification of the target molecules; this helped to improve the flammability, biotoxicity and enrichment effects of the molecule simultaneously. The modified information given by the three-dimension contour maps of the PAEs' flammability and biotoxicity single-effect models indicated the introduction of large groups to DMP was beneficial to the increase of LOI and logLC 50 . So, the DMP derivatives modified by virtue of the comprehensive effect model information could simultaneously meet the requirements of increased flammability and reduced biotoxicity. According to the steric field three-dimension contour map of PAEs' enrichment single-effect model, DMP did not overlap with any color areas, but obviously, the yellow block was bigger. Thus, the introduction of the small group was more conducive to improve the molecular activity; in other words, a large group would reduce the logBCF of the molecules to decrease the enrichment. Therefore, it could be concluded that the modified information of the steric field three-dimension contour maps of the enrichment single-effect model was also consistent with the comprehensive effect.
Polymers 2020, 12, x FOR PEER REVIEW 14 of 21 comprehensive evaluation value. Therefore, in order to reduce the comprehensive evaluation value, larger groups should be introduced into the overlapping part of DMP with three-dimension contour maps during modification of the target molecules; this helped to improve the flammability, biotoxicity and enrichment effects of the molecule simultaneously. The modified information given by the three-dimension contour maps of the PAEs' flammability and biotoxicity single-effect models indicated the introduction of large groups to DMP was beneficial to the increase of LOI and logLC50. So, the DMP derivatives modified by virtue of the comprehensive effect model information could simultaneously meet the requirements of increased flammability and reduced biotoxicity. According to the steric field three-dimension contour map of PAEs' enrichment single-effect model, DMP did not overlap with any color areas, but obviously, the yellow block was bigger. Thus, the introduction of the small group was more conducive to improve the molecular activity; in other words, a large group would reduce the logBCF of the molecules to decrease the enrichment. Therefore, it could be concluded that the modified information of the steric field three-dimension contour maps of the enrichment single-effect model was also consistent with the comprehensive effect. According to the contour maps of the hydrophobic field in Figure 3, the yellow region covered both sides of the DMP molecule in the contour map of PAEs' multi-effect model, which indicated the introduction of the hydrophobic groups was beneficial to reduce the comprehensive evaluation value of the molecule. The contour map of the flammability single-effect model indicated that the introduction of the hydrophobic group at the end of the side chain of the DMP molecule could improve the LOI value, which could enhance the flame retardancy of the molecule. The yellow block was distributed in the side chain of the DMP molecule in the contour map of the biotoxicity model, indicating the introduction of hydrophobic groups could increase the logLC50 of the molecule to decrease the biotoxicity. In the contour map of enrichment model, the volume of the white block was larger, specifying that the introduction of hydrophobic groups was more conducive to improve the molecular activity; this information indicated that the addition of hydrophobic groups during molecular modification could reduce the enrichment.

Mechanism Analysis of PAE Derivatives' Improved Multiple Effects Based on Modified Group Properties
In addition to explaining the molecular modification results from the perspective of three-dimension contour maps, this paper also stated the properties of modified groups by revealing the different improvements on the PAEs' flame retardation, biotoxicity, and enrichment effects; the According to the contour maps of the hydrophobic field in Figure 3, the yellow region covered both sides of the DMP molecule in the contour map of PAEs' multi-effect model, which indicated the introduction of the hydrophobic groups was beneficial to reduce the comprehensive evaluation value of the molecule. The contour map of the flammability single-effect model indicated that the introduction of the hydrophobic group at the end of the side chain of the DMP molecule could improve the LOI value, which could enhance the flame retardancy of the molecule. The yellow block was distributed in the side chain of the DMP molecule in the contour map of the biotoxicity model, indicating the introduction of hydrophobic groups could increase the logLC 50 of the molecule to decrease the biotoxicity. In the contour map of enrichment model, the volume of the white block was larger, specifying that the introduction of hydrophobic groups was more conducive to improve the molecular activity; this information indicated that the addition of hydrophobic groups during molecular modification could reduce the enrichment.

Mechanism Analysis of PAE Derivatives' Improved Multiple Effects Based on Modified Group Properties
In addition to explaining the molecular modification results from the perspective of three-dimension contour maps, this paper also stated the properties of modified groups by revealing the different improvements on the PAEs' flame retardation, biotoxicity, and enrichment effects; the study was performed by introducing groups with different properties. The reliability of the PAEs' multi-effect 3D-QSAR model was proved further.
According to the contour map information of PAEs' flammability, biotoxicity, and enrichment multi-effect 3D-QSAR model, the modified groups were introduced into PAE derivatives having the characteristics of large volume and hydrophobicity. The hydrophobic value of the group was expressed as the logP of the corresponding PAEs derivatives. The larger value of logP indicated higher hydrophobicity. The volume and hydrophobicity values of the modified groups were weighted and coupled to evaluate the strength of the group effect. The weight of the volume and the hydrophobicity of the modified group represented the contribution rate of the steric field and hydrophobic field in PAEs' flammability, biotoxicity, and enrichment multi-effect model, respectively. Further, the coupling values of modified groups' properties were calculated ( Table 9).
The property of each modified group was represented by the average values of the property of coupling values and introduced in Table 9. The average reduction amplitude values of the comprehensive molecular effects, the average increase amplitude values of the flammability effect, the average reduction amplitude values of the biotoxicity, and the average reduction amplitude values of the enrichment of molecules with the same substituted group, were calculated further ( Table 10). As shown in Table 10, the comprehensive effect improvement degree values of these groups (-CH 2 NO 2 , -CH 2 C 6 H 5 , -NO 2 ) were superior to the other introduced groups. The coupling values of each group and the average reduction amplitude values of the comprehensive effect were analyzed by linear regression; in the obtained results, the correlation coefficient R = 0.7516 > 0.7067 (P = 0.05, R = 0.7067), indicated that the modification effect of PAE derivatives was significantly correlated with the volume and hydrophobicity of introduced groups. In addition, the average improvement degree values in flammability, biotoxicity, and enrichment of PAE derivatives were weighted with the weights of 40%: 30%: 30%, and thus the weighted improvement degree values of the comprehensive effect were obtained. The average reduction amplitude values of PAE derivatives' comprehensive effect and the weighted improvement degree values of comprehensive effect were analyzed by linear regression, while the correlation coefficient R of them was 0.8781; in the obtained result, there was a significant correlation at P = 0.01 (R = 0.8343). By virtue of the above results, it was clear that the comprehensive evaluation values of PAEs' flammability, biotoxicity and enrichment effects calculated by the comprehensive evaluation method in this paper could represent the multiple effects of PAEs simultaneously, which confirmed the reliability of the 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment multi-effect.

Conclusions
Using the ideal point method, an improved flammability, biotoxicity, and enrichment multi-effect 3D-QSAR model for PAEs was developed in this paper, and it was successfully applied to the molecular modification of plastic additives PAEs that combined flame retardancy with improvements in biotoxicity and enrichment. Finally, six eco-friendly molecules (DAP-2-CH 2 NO 2 , DAP-1-NO 2 -2-CH 2 C 6 H 5 , DAP-1-NO 2 -2-CH 2 CH 3 , DAP-1-NO 2 -2-CH 2 NO 2 , DAP-1-NO 2 -2-NO 2 , DAP-1-NO 2 -2-OCH 3 , and DAP-2-CH=CH 2 -1-NO 2 ) were obtained, with an increasing rate of LOI value of 0.74-9.12%. Their LOI was between 22-27, which indicates that they were equipped with flame retardancy technically. Their biotoxicity for fish was decreased to a great extent, because the logLC 50 of six molecules was reduced with a range from 127.81-338.22%. Meanwhile, the logBCF values were all less than 3.30 and the molecules would not be accumulated in the organisms. The six molecules' degradability and migration were slightly improved or almost remained unchanged.
In this paper, the reasons for the simultaneous improvement of multiple effects of derivatives were analyzed. Based on the modified information of contour map of PAEs' multi-effect and single-effect 3D-QSAR models for flammability, biotoxicity, and enrichment, it was confirmed that the introduction of large groups and hydrophobic groups may be beneficial to the simultaneous improvement of PAEs' comprehensive effects, and multiple effects of flammability, biotoxicity, and enrichment. In addition, the coupling values of the properties of modified groups were calculated by the volumes of weighting groups and the hydrophobic values; the calculation was performed according to the rate of steric and hydrophobic force field contribution in PAEs' flammability, biotoxicity, and enrichment multi-effect model. The results showed that the coupling values of group properties were directly proportional to the improved degrees of the comprehensive effects of corresponding PAE derivatives, and also correlated significantly (R = 0.7516 > 0.7067 (P = 0.05, R = 0.7067)). The results confirmed the feasibility of the multi-effect 3D-QSAR model of PAEs' flammability, biotoxicity, and enrichment effects.
However, the flame retardancy of PAEs derivatives designed in this paper still had some room to improve, and further attention should be paid to the group types and key quantitative parameters that may affect the flammability performance of PAEs. Funding: This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.