Development of a nano-QSPR model to predict band gaps of spherical metal oxide nanoparticles.

Antibacterial activities and cytotoxicity of metal oxide nanoparticles are determined by their special band structures, which also influence their potential ecological risks. Traditional experimental determination of the band gap is time-consuming, while the accuracy of theoretical computation depends on the selected algorithm, for which higher precision algorithms, being more expensive, can give a more accurate band gap. Therefore, in this study, a quantitative structure-property relationship (QSPR) model, highlighting the influence of crystalline type and material size, was developed to predict the band gap of metal oxide nanoparticles rapidly and accurately. The structural descriptors for metal oxide nanoparticles were generated via quantum chemistry computations, among which heat of formation and beta angle of the unit cell were the most important parameters influencing band gaps. The developed model shows great robustness and predictive ability (R 2 = 0.848, RMSE = 0.378 eV, RMSECV = 0.478 eV, Q EXT 2 = 0.814, RMSEP = 0.408 eV). The current study can assist in screening the ecological risks of existing metal oxide nanoparticles and may act as a reference for newly designed materials. This journal is © The Royal Society of Chemistry.

Introduction

Metal oxide nanoparticles (MONs) are used widely in many fields due to their unique properties.[1] Some of these properties, such as optical and electronic properties, are determined by their special band structure, which can be characterized by three parameters: EC (energy of conduction band), EV (energy of valence band) and Eg (energy of band gap). For example, as a very important optical property, the ability to generate reactive oxygen species (ROS) is tightly related to the band structure of the material. TiO2 and ZnO nanoparticles could generate three types of ROS (·O2−, H2O2,·OH) under UV irradiation while other metal oxides generate only one or two types or even do not generate any type of ROS at all.[2] Different metal oxide nanoparticles also generate different amounts of ROS, and consequently, they have different antibacterial activities or abilities to degrade organic pollutants.[3,4]

Meanwhile, as the wide applications of MONs increase rapidly, their potential ecological risks have drawn more and more attention.[2] One of the main mechanisms of oxidative damage caused by MONs is the overlap of conduction band energy (EC) levels with the cellular redox potential. The overlap of EC with the redox potential leads to the transfer of electrons from reductive substances, which then generate reactive oxygen species (ROS) which could eventually damage cells.[3-5] Therefore, the EC value is a meaningful parameter for the assessment of the ecological risks of MONs. EC can be calculated by means of eqn (1): EC refers to the conduction band energy; Eg is the band gap energy and χoxide is the absolute electronegativity of the metal oxide.[6]

The absolute electronegativity can be calculated by a set of equations reported by Portier et al.[7] Values of the band gaps of MONs can be obtained from experimental determination and from estimation.

The experimental determination of band gaps requires diffuse reflectance UV-vis spectroscopic analysis, which can accurately measure the band gaps of MONs. However, considering the increase of emerging MONs both in amount and types as well as taking into account that band gaps of MONs vary when crystalline types, sizes and shapes are different,[8] experimental determinations of band gaps are inefficient and cannot provide any guide at all for the design of novel MONs that are safe with regard to their ecological impacts.

Estimation methods were proposed, including theoretical computation and empirical equation to obtain band gaps efficiently and conveniently. As for theoretical computation, its inaccuracy is the main obstacle which undermines its authorization and application domain. For example, the band gaps of solids could be underestimated with an accuracy of 30–100% in a typical density functional theory calculation.[9] Building a cluster of unit cells was tried,[10] but computation for a cluster commensurate with real material requires vast investment of computational resources. As for empirical estimation, Portier et al.,[11] demonstrated the correlation between band gap and thermodynamic properties, and put forward an empirical equation to estimate band gap of bulk metal oxides by using enthalpies of formation (eqn (2)).

In the eqn (3), refers to the standard enthalpy of formation of the oxide in cal mol−1 and the coefficient A is determined by cation ion of the oxide, N is Avogadro's number, ne refers to the number of valence electron. For s-block, p-block, d-block and f-block of metal elements, the reported A values are 0.80, 1.35, 1.00 and 0.80 respectively.[11] Using this empirical equation, the band gaps of binary oxides can be calculated very quickly. But this equation ignores the effects of size differences. Therefore, it is required to develop a prediction method to obtain the band gap of MONs with different crystalline type and particle size accurately and efficiently.

Quantitative structure–property relationships (QSPRs) have been used to predict various properties of chemicals.[12,13] QSPRs can be applied fast, they are low-cost and efficient compared with other methods. For MONs, several QSPRs models were reported to predict photo-induced toxicity of MONs, including toxicity of MONs to E. coli and to human keratinocyte cells as well as the cell uptake rate of MONs.[14-17] Wyrzykowska et al. developed a nano-QSPR model to predict the zeta potential of MONs in KCl aqueous solution of MONs.[18]

In this study, we developed a nano-QSPR model, considering different crystalline types and sizes, to predict band gaps of MONs. Besides, band gaps were also calculated by two other estimation methods to compare the predictive ability of different methods. From a wider perspective, this study could help screen the potential ecological risk of MONs and provide reference for newly designed MONs.

Methods

Data set for band gaps

Experimental band gaps of MONs were obtained by collection from literature[3,19-50] and by experimental determination.

To control variables that influence band gaps of MONs as much as possible and to ensure the consistency of collected band gaps, a qualified data point from literature must meet the following demands: size of material is in the nanometer scale; particles must be spherical or approximate a spherical shape; particles have a single chemical makeup without any chemical modification on the surface; characterizations including X-ray diffraction, which verifies the crystalline structure of MONs, and UV-vis spectroscopic analysis are required.

Seven kinds of MONs (5 nm anatase-TiO2, 25 nm anatase-TiO2, 40 nm anatase-TiO2, 100 nm anatase-TiO2, 40 nm rutile-TiO2, 100 nm rutile-TiO2 and 50 nm ZnO) bought from Aladdin (http://www.aladdin-e.com/) were prepared for UV-vis spectroscopic analysis using UV-550 spectrophotometer produced by JASCO Corporation. To minimize possible errors, each determination for each material was repeated three times and the eventual outcome is the arithmetic mean of parallel experiments. The raw data of UV-vis spectroscopic analysis were transformed according to the Kubelka–Munk equation.[51]

The total data set of band gaps of MONs was split randomly into a training set and a validation set with a ratio of 3 : 1 (30 : 10).

Date pre-processing

Previous studies indicated that the band gap of relatively large particles (size > 10 nm) is not affected by size whereas the quantum confinement effect only occurs when the size of nanoparticle is smaller than approximately 10 nm, consequently leading to a larger value of band gap.[52,53] Therefore, the value of the band gaps of materials with size smaller than 10 nm was used directly. In case of material with a size larger than 10 nm, the average size and the average band gap were calculated and included in either the training set or the validation set.

Descriptors

The descriptors used to develop the nano-QSPR were derived from theoretical calculation, as well as from literature and experimental characterization. (1) Theoretical descriptors were calculated by using the Vienna Ab initio Simulation Package (VASP) 5.4.1 (ref. 54) and MOPAC2016 (ref. 55) (2) reported descriptors of MONs in literatures, such as the ratio of metal atoms to oxygen atoms and weight of atoms, were included in this study, remaining to be selected further.[4,56] (3) Experimental characterization outcomes such as size, crystalline type information and the number of unit cells contained by the material were included.

To consider the influence of crystalline types to band gaps of materials in calculation, input files should specify crystalline type, which must be consistent with corresponding band gap of every material. The XRD outcomes or joint committee powder diffraction standards (JCPDS) number was collected from the same literature from which the value of the band gap of MONs was collected. MDI JADE 5.0 (ref. 57) was utilized to analyze collected XRD patterns and to match the corresponding standard XRD cards, which could specify crystalline type. Based on the crystalline type information derived from standard cards, the corresponding crystalline file in common intermediate format (CIF) was obtained in FindIt[58] before being transformed by Open Babel 2.4.1 (ref. 59) into standard input files. Single cell units were optimized by the general gradient approximate (GGA) of VASP before descriptor computation and the error of optimization for structure was verified to be within 5%. Descriptors of VASP computation were obtained from POTCAR and OUTCAR after structure optimization. The PM7 method was used to generate various descriptors in MOPAC computation. In order to validate the acceptability of descriptors derived from the computation being based on single cell units, important descriptors selected in model were compared with their experimental values.

In order to ensure that the correlation between descriptors and band gaps is linear, several parameters were processed before they were suited as descriptor. Heat of formation of the unit cell, computed by MOPAC, was transformed according to a previously reported empirical equation.[11] The reciprocal value of the squared diameter was processed as a descriptor according to the Brus equation (eqn (4)).[52]

QSPR modeling

Partial least square (PLS) regression was utilized to develop the QSPR model. The SIMCA 13.0 software package was used to provide PLS as well as the subsequent analysis.[60] The PLS approach was established around 1975 by Herman Wold for the modelling of complicated data sets.[61] The basic principle of PLS could be interpreted as PLS being the combination of principal components analysis (PCA), correlation analysis (CA) and multiple linear regression (MLR).[62] Unlike MLR, PLS can analyze data with strongly collinear, noisy and numerous X variables.[63] In our study, the amount of descriptors is relatively large and the correlation between various descriptors is difficult to demonstrate. Therefore, PLS was used to develop this model.

A genetic algorithm (GA), developed by means of the R program,[64] was used to select descriptors for developing prediction models. The number of descriptors in the model was further decreased by using a method reported in literature:[65] Based on the values of the variable influence for the projection (VIP), descriptors with lower VIP values were deleted from the model until the model shows best validation performance.

Model characterization

In order to ensure the reliability and predictive ability of the developed model, various analyses were done to characterize the performance of the model. The correlation coefficient R2 between predicted and observed Eg and the root mean standard error of estimation (RMSEE) measure the fitting performance. The performance of internal validation could be indicated by the root mean standard error of cross validation (RMSECV), which is based on the leave-one-out (LOO) test. As for external validation, the predictive power of the developed QSPR model was estimated by means of calculating the root mean standard error of prediction (RMSEP) and QEXT2 defined by means of eqn (5), where y and ŷ are the measured and predicted value of Eg in the validation set, respectively; ȳ is the average value of Eg for the training set. A value of QEXT2 larger than 0.5 indicates a good predictive capability of the model and a value of QEXT2 larger than 0.9 indicates excellent predictive ability.

According to the OECD guideline,[66] the applicability domain of the developed model was determined by the method reported by K. Roy et al.[67] Firstly, descriptors in training set were normalized. If all of the normalized descriptors for a material are larger than 3, that material can be judged as an outlier. On the contrary, a material is treated as non-outlier if all of its normalized descriptors are smaller than 3. For a material who has some descriptors larger than 3 but the others not, it will be further investigated.

Results and discussion

Data set of band gaps

Through data collection, 91 values of band gaps of 22 different MONs were collected.[3,19-50] Band gaps ranged from 1.85 eV to 6.07 eV with an average of 3.60 eV and a standard deviation of 0.8 eV. The particles diameters ranged from 2.5 nm to 100 nm with an average of 26.0 nm and a standard deviation of 24.8 nm. All collected data are shown in Table S1† and the data set processed according to particles size for model development is displayed in Table 1.

Metal oxide nanoparticles and their structural descriptorsa

No.	Metal oxide	a	b	c	α	β	γ	Size (nm)	E exp. (eV)	E pred. (eV)	Set
1	Al2O3	4.81	4.81	13.12	90	90	120	76.0	5.97	5.74	T
2	CeO2	5.46	5.46	5.46	90	90	90	12.0	3.39	3.74	T
3								3.6	3.68	3.82	T
4								2.6	3.44	3.89	T
5	Cr2O3	4.58	4.58	14.72	90	90	120	26.5	4.22	4.11	T
6								38.2	3.70	4.11	V
7	Cu2O	4.31	4.31	4.31	90	90	90	4.0	2.90	2.35	T
8								9.0	2.50	2.29	V
9								6.0	2.83	2.31	V
10								20.0	2.04	2.28	V
11	Fe2O3	4.74	4.74	13.49	90	90	120	100.0	2.38	2.96	V
12	Fe3O4	8.04	8.04	8.04	90	90	90	12.0	1.85	2.61	T
13	Ga2O3	12.50	3.10	5.92	90	103.7	90	65.0	4.85	4.81	T
14								2.5	5.10	5.00	T
15	HfO2	5.14	5.19	5.31	90	99.2	90	17.0	6.07	5.66	T
16	In2O3	10.35	10.35	10.35	90	90	90	17.0	3.85	4.37	T
17	La2O3	3.94	3.94	6.18	90	90	120	30.0	4.98	4.30	T
18	MgO	4.25	4.25	4.25	90	90	90	7.0	4.27	4.54	T
19	Mn2O3	9.03	9.03	9.03	90	90	90	30.2	3.27	3.38	T
20	NiO	4.16	4.16	4.16	90	90	90	3.5	3.67	3.62	T
21								4.6	3.63	3.58	T
22								5.5	3.62	3.57	V
23								10.1	3.62	3.54	T
24								12.4	3.61	3.54	T
25								22.0	3.27	3.53	T
26	Sb2O3	5.18	16.61	5.51	90	90	90	11.8	4.49	4.77	T
27	SnO2	4.83	4.83	3.24	90	90	90	4.5	4.20	3.80	T
28								4.0	4.10	3.81	V
29								5.0	4.21	3.79	T
30								5.2	4.20	3.79	V
31								3.7	4.33	3.82	T
32								30.0	3.74	3.75	T
33	TiO2-a	3.82	3.82	9.70	90	90	90	17.3	3.25	3.44	T
34	TiO2-b	9.28	5.52	5.19	90	90	90	7.9	3.48	3.85	T
35								18.9	3.11	3.83	V
36	TiO2-r	4.66	4.66	2.97	90	90	90	70.0	3.00	3.00	T
37	WO3	5.39	5.36	7.84	90	91.8	90	42.0	2.77	2.98	T
38	Y2O3	10.65	10.65	10.65	90	90	90	14.6	5.30	5.38	V
39	ZnO	3.29	3.29	5.31	90	90	120	35.0	3.24	2.62	T
40	ZrO2	5.19	5.24	5.38	90	99.2	90	40.0	5.04	5.36	T

T: training set; V: validation set; Eexp.: experimental gaps; Epred.: predicted gaps.

The distribution of band gaps values for MON with multiple sizes was investigated in more detail. First of all, the arithmetic means and standard deviations of the band gaps were calculated: CeO2 (3.50 ± 0.13 eV), Cr2O3 (3.82 ± 0.24 eV), Cu2O (2.57 ± 0.34 eV), Fe2O3 (2.38 ± 0.33 eV), Ga2O3 (4.93 ± 0.12 eV), NiO (3.43 ± 0.19 eV), SnO2 (3.93 ± 0.38 eV), TiO2 (3.25 ± 0.03 eV), ZnO (3.25 ± 0.08 eV). Fig. 1a and b show the distribution of band gaps of NiO and SnO2 respectively as a function of their squared diameter. It can be observed clearly from these plot that band gaps increase as size decreases to values below 10 nm, indicating the occurrence of a quantum confinement effect around 10 nm. However, ambiguity does exist. For example, band gaps increase was not observed for some particles smaller than 10 nm. Besides, large variance of values of band gaps were observed, especially for SnO2 larger than 10 nm. NiO with size of 10.1 nm and 12.4 nm were treated as two single points although their sizes are larger than 10 nm because they display significant difference from other data with size larger than 10 nm. It is the same with Cr2O3. These phenomena might due to these materials are from different literatures, in which materials were synthesized by different methods in different conditions. These differences are hard to characterize because of limited information. Therefore, using the average value of data points with size larger than 10 nm could avoid this unexplainable difference and will minimize the uncertainty of data generated in different laboratories.

Fig. 1

(a) Band gaps of NiO versus the squared diameter of the particles. (b) Band gaps of SnO2versus the squared diameter of the particles. The dotted line marks the value of 100 nm2, which means that data points in the zone left to the dotted line refers to particles with a size smaller than 10 nm and the other data points relate to particles that are larger than 10 nm.

Model development and evaluation

Through PLS regression and descriptor selection, the model was obtained as given in eqn (6) (standard error of regression coefficients are listed in Table S5†). The coefficient of every descriptor in the reported equation is unscaled in order to be convenient for subsequent use. Fig. 2 was plotted to visualize the observed experimental Eg values vs. predicted values of Eg as calculated by eqn (6). The correlation coefficient R2 and the RMSEE of the training set containing 30 values of band gaps of MONs are 0.848 and 0.378 eV, indicating an excellent goodness-of-fit. QEXT2 and RMSEP are 0.814 and 0.408 eV respectively, indicating good predictive ability. The developed model contains three principal components (N : A = 10 : 1, N: number of data in training set, A: number of principal components), which accounts for 85 percent explanation of the total variance. The physical meaning and the loading weights of the descriptors in the first two principal components are displayed in Table 2. In the first principal component, heat of formation, beta angle of unit cell and the Thomas Fermi vector outweigh the other descriptors. For the second principle component, Fermi energy, heat of formation and the ratio of metal atoms to oxygen atoms are the most significant descriptors. The VIP values of selected descriptors are plotted in Fig. 3. The VIPs of heat of formation, beta angle of unit cell and Thomas Fermi vector are larger than 1, indicating that these descriptors explain most of the variance of band gaps. Computed heat of formation was compared with the experimental ones in Table S4 and Fig. S1† to validate the acceptability of computed descriptors derived from small molecular model computation, which displayed acceptable fitting performance.N = 30, A = 3, R2 = 0.848, RMSEE = 0.378 eV, RMSECV = 0.478 eV, QEXT2 = 0.814, RMSEP = 0.408 eV.

Fig. 2

Observed vs. predicted Eg values.

Meaning and loading weight of descriptors in the QSPR modela

Descriptor	Physical meaning	Source	w*c [1]	w*c [2]
H F	Heat of formation	MOPAC	0.65	0.57
BETA	Beta angle of unit cell	—	0.51	0.34
D −2	Reciprocal of square diameter	—	0.04	0.03
V 2	Length of second vector of unit cell	VASP	−0.06	−0.12
E Fermi	Energy of Fermi level	VASP	0.03	0.72
TFW	Thomas Fermi vector	VASP	−0.35	0.14
R	Ratio of metal atoms to oxygen atoms	—	−0.08	0.45
E T	Total energy	MOPAC	0.26	0.12
DENC	Electron–electron repulsion energy	VASP	0.28	−0.04
XCENC	Exchange and correlation energy	VASP	−0.20	0.14

w*c [1]: loading weight of component [1] w*c [2]: loading weight of component [2].

Fig. 3

Variable influence for the projection (VIP) plot of selected descriptors.

The applicability domain of the developed model was determined by the reported method mentioned above. Data points from both training set and validation set were investigated. The detailed information of the calculation process is provided in Table S2.† There are 5 and 4 data points with at least one outlier of descriptor in the training set and the validation set respectively. But there is only one material (Y2O3 in validation set) that was judged as an outlier, mainly because its heat of formation deviates far from that of most metal oxides in this model.

Mechanistic interpretation

Heat of formation and beta angle of unit cell are the most significant descriptors. The correlation between band gap and heat of formation has already been demonstrated by previous studies, that is: bulk metal oxides with larger heat of formation tend to have a larger band gap.[68-70] In this study, this correlation for metal oxide nanoparticles was also observed (p < 0.001), which means that this relation remains irrespective of a decrease of material size. It needs to be notified that in this study, heat of formation was calculated by MOPAC PM7. This calculation may not be as accurate as experimental determination but provides a quick method for obtaining heat of formation. Crystalline types could be reflected by not only direct variables characterizing crystalline type such as beta angle and length of vector but also by indirect variables like the E-Fermi and the Thomas Fermi vector because these parameters were calculated based on specific crystalline types. The presence of the beta angle (p < 0.05) and the length of unit cell in the developed model indicates that crystalline type is an important parameter influencing the band gap of MONs. Beta angle and length of vector contribute to the interference of the wave vector of electrons and consequently determine the density of electrons and eventually influence band gaps.

The VIP of materials size is not large, which coincides with the Brus equation and the computed standard deviation of materials with multiple sizes. Although the effect of size has been taken into consideration during data pre-processing, it needs to be noted that some variables which may influence the band gap of MONs, were not contained in the developed model. For example, the annealing temperature, time and the ratio of different precursor materials might also affect the value of the band gap. However, these variables are difficult to describe because of the numerous types of MONs and numerous methods to synthesize nanoparticles. Different experimental conditions may lead to unknown differences in the composition of materials, which need to be characterized further.

Model comparison

The developed QSPR model was compared with empirical equation. Band gaps were calculated with the computed and experimental heat of formation according to the empirical equations.[11] The results of comparison are shown in Table 3.

Predictive performance of the developed model vs. empirical equation

Predictive method	R 2	RMSE (eV)	Reference
QSPR model	0.848	0.378	This study
Empirical equation based on computed heat of formation	0.490	2.24	Portier et al., 2001
Empirical equation based on experimental heat of formation	0.571	1.58	Portier et al., 2001

The empirical equation using experimental heat of formation to predict band gap, shows better performance than the empirical equation which used computed heat of formation. However, R2 did not change significantly while there was a decrease in RMSE. During characterization the performance of these two empirical methods, it is found that these two methods underestimate or overestimate band gaps. Our QSPR model developed by PLS performs better than previous empirical equations.

Conclusions

A nano-QSPR model, covering various crystalline types and material sizes, was developed to predict band gaps for metal oxide nanoparticles. The established QSPR model performs well in terms of goodness-of-fit, robustness and predictive ability. Heat of formation and the beta angle of unit cell are the most important parameters influencing the band gaps of metal oxide nanoparticles. Compared with empirical equation related to the heat of formation, the developed QSPR model had a larger value of R2 and a smaller RMSE. This newly developed QSPR model could thus contribute significantly to the rapid screening of the ecological risk of existing MONs and could assist the design of novel MONs that are environmentally safe.

Conflicts of interest

The author declares no competing financial interest.