research papers
The generalized F constraint in the maximumentropy method – a study on simulated data
^{a}Laboratory of Crystallography, University of Bayreuth, Germany
^{*}Correspondence email: palat@unibayreuth.de
One of the classical problems in the application of the maximumentropy method (MEM) to electrondensity reconstructions is the uneven distribution of the normalized residuals of the structure factors of the resulting electron density. This distribution does not correspond to the expected Gaussian distribution and it leads to erroneous features in the MEM reconstructions. It is shown that the classical constraint is only one of many possible constraints, and that it is too weak to restrict the resulting distribution to the expected Gaussian shape. It is proposed that constraints should be used that are based on the higherorder central moments of the distribution of the structurefactor residuals. In this work, the influence of different constraints on the quality of the MEM reconstruction is investigated. It is proposed that the use of a combined constraint on more than one central moment simultaneously would lead to again improved results. Oxalic acid dihydrate was chosen as model structure, from which several data sets with different resolutions and different levels of noise were calculated and subsequently used in the MEM. The results clearly show that the use of different constraints leads to significantly improved results.
Keywords: maximumentropy method; electron density; oxalic acid.
1. Introduction
The maximumentropy method (MEM) is used as a powerful tool for a modelfree image reconstruction in many scientific applications (von der Linden et al., 1998). In crystallography, one particular application is the investigation of the electron density in the crystal structure. After the first promising applications in this field (Collins, 1982; Sakata & Sato, 1990), several warnings concerning the reliability and possible pathologies of the method appeared (Jauch, 1994; de Vries et al., 1996). One of the obvious problems was that the distribution of the normalized residuals of the structure factors strongly deviated from the expected Gaussian distribution. Some of the reflections – usually strong reflections at low angles – had very large , while the others were fitted almost exactly. The large deviation of the histogram of from the Gaussian distribution was responsible for unphysical features in the corresponding electron density. A solution to this problem was proposed by de Vries et al. (1994), who employed an ad hoc weighting scheme within the classical constraint. However, a theoretical basis for this weighting scheme does not exist.
Here we propose new constraints based on the higherorder central moments of the distribution of . We show that the use of these constraints produces results with better distributions of and with less artifacts in the reconstructed electron density than the classical constraint.
The method is tested against data sets of various resolutions and with various noise levels that were computed for a known electron density of oxalic acid dihydrate.
2. The method
The basic principle of the MEM is that the optimal image is defined to be the image with the maximum value of the S, while one or more constraints of the type are fulfilled. For our purposes, the image is the electron density () in the which is defined by its values on a grid of N_{p} = N_{1}×N_{2}×N_{3} points. The is defined as
functionalwhere the values define the prior or reference electron density. For an overview of the crystallographic applications of the MEM, see Gilmore (1996). The constraints should be selected so as to define which image is in agreement with the observed data. The first reasonable constraint is the normalization of to the expected number of electrons per unitcell volume:
Traditionally, the constraint to the scattering data is the leastsquares likelihood criterion , with
where the summation runs over all measured structure factors N_{F}. This definition of the constraint cannot be used directly, since it does not contain the information about the phases of the structure factors and does not lead to convergence. Therefore, the socalled F constraint is usually employed:
The value of C_{F} depends on both the amplitudes and phases of F_{obs} and F_{MEM}. C_{F} is minimal if the phases of all F_{obs}^{i} are equal to the corresponding F_{MEM}^{i}. In that case, .
The use of the statistics (and its phased modification in the C_{F} constraint) is based on an assumption that the experimental errors on are random with a Gaussian distribution:
where is a sample of the random variable with normalized Gaussian distribution. Since the resulting electron density should be the best estimate of the true density, the corresponding calculated structure factors F_{MEM} should be the best estimate of F_{true} and the distribution of the normalized residuals should be Gaussian too.
It is obvious that the Gaussian distribution of errors does imply the validity of the (or C_{F}) constraint, but not vice versa. Constraining only is not sufficient to ensure the proper Gaussian form of the resulting error distribution.
A probability distribution of a random variable x is characterized by the values of its central moments m_{n}. For the normalized Gaussian distribution, the central moments are
The values of the moments of odd order are all zero and the moments of even order are:
In the case of N samples of the variable x, the central moments m_{n} can be computed from
It follows from (3) and (8) that is the m_{2} central moment of the distribution of . Thus, the concept of generalized F constraint can be introduced, with F_{2} referring to the classical constraint on the secondorder moment, and with F_{n} defining a constraint based on the moment of order n:
Only the constraints with n even restrict the width of the histogram, constraints with n odd are sensitive only to the symmetry of the distribution with respect to the origin. Therefore, only the constraints with n even are used in this work.
It has been suggested that more simultaneous constraints (up to the number of independent observations) of the form could be used instead of the single constraint (Carvalho et al., 1996). This requires some additional criterion for defining the point of convergence and strongly restricts the role of the MEM as the noise filter. We suggest that the use of several F_{n} constraints simultaneously is the proper way to handle noisy data, since the expected shape of the histogram is the only information about the noise that is available. However, the available algorithms do not allow such a generalization. Therefore, in the present stage of the work, the influence of different choices of a single constraint based on (9) on the result of MEM was investigated.
3. Computational details
The method was tested on the structure of oxalic acid dihydrate. The main reason for this choice was that this compound has become a kind of standard for chargedensity studies. In addition, the structure of oxalic acid dihydrate is very suitable for this type of work, since it is centrosymmetric and the central molecule is planar. This allows an easy interpretation of the majority of the features using only one section of the electron density. The basic characteristics of the structure are summarized in Table 1.

First, the electron density of the procrystal structure (superposition of independent atoms, ) was created. This was done by a method due to Papoular et al. (2002). The analytical approximation to spherical atomic scattering factors (Su & Coppens, 1997) for each atom of the structure was multiplied by the anisotropic displacement factor of that atom. The resulting threedimensional distribution in was then transformed by means of the analytical Fourier transform to obtain the electron density of that atom. The density was sampled on the 64 ×32 ×128 pixel grid, which corresponds to a pixel size of approximately 0.1×0.1 ×0.1 Å. The positional and displacement parameters from the due to Šlouf (2001) were used. The electron densities of the individual atoms were then summed to obtain . The `true' electron density was then constructed by summing with the dynamic deformation density , as determined by the multipole of Šlouf (2001) (Fig. 1a). This caused 1.65% of the pixels of the resulting electron density to be negative. The lowest density was −0.021 e Å^{−3}. The negative areas were located in the lowdensity intermolecular regions. This unphysical feature probably originates from the inaccuracy of the multipole expansion in these very low density regions. The MEM cannot handle these negative regions and very low density regions increase the dynamic ratio of the electron density inadequately. Therefore, the pixels with e Å^{−3} were set to 0.005 e Å^{−3}. 2.45% of the pixels were corrected.
The electron density obtained by this procedure is certainly not the true electron density of oxalic acid dihydrate. The analytical approximation used in the first step is not absolutely accurate and the structure parameters and multipole deformation density can contain a substantial degree of inaccuracy, too. However, this model of electron density is good enough to be used as the reference electron density for MaxEnt calculations and will be denoted as (Fig. 1b).
The structure factors corresponding to the original map were calculated by means of a numerical Fourier transform. To investigate the influence of noise and resolution on the quality of the MEM reconstruction, 16 different data sets were created. The value is used as a measure of resolution in this paper. It was chosen to be 0.5, 0.75, 1.0 and 1.25 Å^{−1} for the respective data sets, and for each resolution four different levels of Gaussian noise were added to the calculated structure factors. To simulate the error distribution in real experimental data, were calculated from
where defines the noise level, simulates the influence of nonzero background and p is the commonly used instability factor. The noisy `observed' structure factors were then calculated to fulfil the equation
Here, is a random variable with normalized Gaussian probability distribution. Three different nonzero noise levels were created this way. The noiseless data sets at each resolution were included for checking purposes. Although the structure factors in the noiseless data sets were exact, which means they should be assigned a zero standard deviation, this is not possible owing to the nature of the constraints [equation (9)]. Therefore, the value of was set to 0.005 for all structure factors so as to be low enough and to allow the computations to finish in a reasonable time. The parameters of different noise levels and resolutions are summarized in Table 2 and Fig. 2.

It is interesting to compare the phases of structure factors corresponding to with the phases corresponding to . In the present case, which is representative for investigations of accurate electron densities, the amount of the unknown structure is minute and the phases of the true structure factors are very well estimated by the phases of the structure factors of . Among all 4029 structure factors, up to Å^{−1}, only nine have different phases for and . Moreover, equation (11) allows for changes of phases between F_{obs} and F_{true}. As a consequence of the introduction of the noise, there have been many more phases changed in each noisy data set than nine. Thus, the results presented here are not influenced by the preliminary multipole and can be regarded as being obtained using just the standard refinement.
We have developed our own computer program BayMEM (first version by Schneider, 2001) for the application of the MEM in chargedensity analysis. This program is designed to work in general ndimensional space to allow computations of the MEM electron density of incommensurately modulated structures, but can be used for standard threedimensional structures too without any restrictions. BayMEM can use both the algorithm of Sakata & Sato (1990) and the MEMSys5 package (Gull & Skilling, 1999). The program was extended to deal with the generalized F constraint. For the present study, the algorithm by Sakata & Sato (1990) was used.
The following characteristics are used to compare the quality of the MaxEnt reconstructions:

The quality of the MEM reconstructions can be compared with Fourier maps. The Fourier transform of the observed structure factors with calculated phases results in an electron density () that can be compared with as obtained with the uniform prior. Inspection of shows that the noise is much larger than in . This is quantified by the C values (Table 4).
The classical method to derive information about electron densities beyond the model is the difference Fourier. We have computed the difference Fourier for F_{obs}F_{pro} (). To be able to compare the result with , we have added to . Again, the noise in is significantly larger than in (Table 5).
4. Results and discussion
4.1. The uniform prior
In the first series of calculations, a uniform electron density was used as prior. The dominating structure of is the oscillatory electron density around each atomic position (Fig. 3). Its presence is independent of the constraint and of the noise level. However, at high noise levels these features are partly camouflaged by the noise of itself. The oscillations are most pronounced at the zero noise level. Clearly, this effect is a demonstration of the seriestermination error intrinsically present in the method, as pointed out already by Jauch (1994) and later discussed in detail by Roversi et al. (1998). The present results show the extent of this effect and its dependence on the resolution of the data set. The amplitude of the artifacts decreases with resolution, but even at resolution 1.25 Å remains significant (Fig. 3, Table 3). Further lowering of the artifacts by increasing the resolution is in practice not possible due to the experimental limitations. Possible ways to overcome this problem are summarized in §5.

The obtained for different noise levels and different resolutions is characterized by the C values (Table 4), by the shapes of the histograms of (Fig. 4), and by the values of the central moments of the distribution of (Fig. 5). The following conclusions can be made based upon the table and the figures:
). The waviness of the lowdensity contours in is suppressed, the overall amount of the residual structure in decreases. It should be noted that the total density maps do not give sufficient insight into the accuracy of the result and cannot be used as a single criterion of the quality of the MaxEnt reconstruction. This can be seen from the comparison of the total and difference maps (Fig. 6). The largest errors occur in the medium and high density levels, where the total density map seems to be smooth and well behaved. This is especially true for the lowresolution maps, which seem to be smooth at first sight, but which exhibit large differences in comparison to the original map.Despite the significant improvement of the MEM reconstructions obtained with the constraints on the higherorder moments, the quality of the reconstructions using the static weighting was in our case even better (Table 4). This surprising effectiveness of the idea of the static weighting suggests that there might exist some fundamental reason for it. A closer investigation of possible theoretical foundations of this type of weighting is desirable.
The systematic investigation of the large number of different data sets allows one to make some general conclusions about the influence of the noise and the resolution on the quality of the result. The expected improvement of the C factors with decreasing noise level is clearly visible. The improvement with the increasing resolution is visible, too, but not as an absolute rule (compare C values of n3r1.00 and n3r1.25, n2r1.00 and n2r1.25 in Table 4). This can be correlated with Fig. 2. The larger the fraction of unobserved reflections present in the outer shell, the smaller is the amount of information it contains. In the data sets with the high noise level, almost all reflections in the outer shells are lessthan's, and they cannot contribute to the improvement of the MEM reconstruction.
4.2. The procrystal prior
In the second series of calculations, the procrystal electron density was used as prior. The summary of the resulting C values is given in Table 5. The deformation density obtained with data sets n2r1.00 and n1r0.75 is shown in Fig. 7. We believe that these examples are quite close to the data sets obtainable in practice.

As expected, the artifacts are strongly reduced and visible only in the vicinity of the atomic center. The deformation density resembles the true deformation density quite well even for the medium noise level. The differences in C factors among the different F_{n} constraints and the different static weighting are much smaller than in the case of the uniform prior, but they are still significant, especially for the low noise levels.
With increasing noise level, the outer shells of structure factors contain so much noise that it masks their statistical difference from the prior structure factors. Such reflections do not improve the result and can even lead to a slightly worse (compare Table 5 and Fig. 2). In an extreme case – noise level 3 – the reflections do not provide any additional information at all and is almost identical with the prior. In other words, the MEM indicates that the data do not contain any evidence for deviation from the prior.
The results confirm that, with procrystal prior information, the MEM is able to reveal the deformation electron density even from the mediumresolution data, provided they are sufficiently accurate.
5. Conclusions
The intrinsic presence of the seriestermination effect in the crystallographic applications of the MEM is demonstrated. The extent of this effect depends on the resolution of the data set and on the kind of prior electron density. For the uniform prior, the artifacts are significantly higher than the bonding electrondensity level and make this version of the MEM unsuitable for investigation of fine features in the electron density. Nevertheless, it is still a useful method for investigation of more robust features like anharmonic atomic movement or disorder (Bagautdinov et al., 1998; Dinnebier et al., 1999; Wang et al., 2001).
The procrystal prior electron density lowers the artifacts and the reconstructions with this prior contain the information about the fine features of the electron density. Further lowering of the artifacts could probably be achieved with the twochannel MEM (Papoular et al., 1996) or with the valenceonly MEM proposed by Roversi et al. (1998). The latter method uses the refined structure parameters to create a core electrondensity fragment, which is then considered to be known and is not included in the MaxEnt optimization. Only the valence electron density is modified. However, this method is of practical use only for extremely accurate data from simple structures, since it relies on the knowledge of the temperature parameters, which are often inaccurate and correlated with systematic errors in the data sets.
The use of the generalized F constraint dramatically improves the quality of the MEM results. The selection criterion for the proper order is the best coincidence of the histogram with the expected Gaussian distribution. From our experience, the order 4 or 6 gives the best result.
Static weighting still gives better results than the nonweighted F_{n} constraints. But this type of weighting lacks any theoretical foundation, and the choice of the best weighting is very data set dependent (Yamamoto et al., 1996). On the other hand, the constraints based on the expected moments of the distribution of have a clear interpretation. One can expect that the new algorithms that will allow the simultaneous use of several constraints in the MEM will again lead to improved results.
One more advantage of the higherorder F constraints in comparison to the classical F_{2} constraint or static weighting is faster convergence, which makes the computation time significantly shorter.
Acknowledgements
Financial support by the Deutsche Forschungsgemeinschaft is gratefully acknowledged.
References
Bagautdinov, B., Luedecke, J., Schneider, M. & van Smaalen, S. (1998). Acta Cryst. B54, 626–634. Web of Science CrossRef CAS IUCr Journals
Carvalho, C., Hashizume, H., Stevenson, A. & Robinson, I. (1996). Physica (Utrecht), B221, 469–486. CrossRef Web of Science
Collins, D. M. (1982). Nature (London), 298, 49–51. CrossRef CAS Web of Science
Dinnebier, R. E., Schneider, M., van Smaalen, S., Olbrich, F. & Behrens, U. (1999). Acta Cryst. B55, 35–44. Web of Science CSD CrossRef CAS IUCr Journals
Gilmore, C. J. (1996). Acta Cryst. A52, 561–589. CrossRef CAS Web of Science IUCr Journals
Gull, S. F. & Skilling, J. (1999). MemSys5 v1.2 Program Package. Suffolk, United Kingdom.
Jauch, W. (1994). Acta Cryst. A50, 650–652. CrossRef CAS Web of Science IUCr Journals
Linden, W. von der, Dose, V., Fisher, R. & Preuss, R. (1998). Editors. Maximum Entropy & Bayesian Methods. Dordrecht: Kluwer Academic Publishers.
Papoular, R. J., Vekhter, Y. & Coppens, P. (1996). Acta Cryst. A52, 397–407. CrossRef CAS Web of Science IUCr Journals
Papoular, R. J., Collin, G., Colson, D. & Viallet, V. (2002). In Proccedings of the 21st Workshop on Bayesian Inference and Maximum Entropy Methods in Science and Engineering, edited by R. Fry. Melville, NY: American Institute of Physics. To be published.
Roversi, P., Irwin, J. J. & Bricogne, G. (1998). Acta Cryst. A54, 971–996. Web of Science CrossRef CAS IUCr Journals
Sakata, M. & Sato, M. (1990). Acta Cryst. A46, 263–270. CrossRef CAS Web of Science IUCr Journals
Schneider, M. (2001). PhD thesis, University of Bayreuth, Germany.
Šlouf, M. (2001). PhD Thesis, Charles University, Prague, Czech Republic.
Su, Z. & Coppens, P. (1997). Acta Cryst. A53, 749–762. CrossRef CAS Web of Science IUCr Journals
Vries, R. Y. de, Briels, W. J. & Feil, D. (1994). Acta Cryst. A50, 383–391. CrossRef Web of Science IUCr Journals
Vries, R. Y. de, Briels, W. J. & Feil, D. (1996). Phys. Rev. Let. 77, 1719–1722. CrossRef Web of Science
Wang, C.R., Jai, T., Tomiyama, T., Yoshida, T., Kobayashi, Y., Nishibori, E., Takata, M., Sakata, M. & Shinohara, H. (2001). Angew. Chem. Int. Ed. Engl. 40/2, 397–399. Web of Science CrossRef
Yamamoto, K., Takahashi, Y., Ohshima, K., Okamura, F. P. & Yukino, K. (1996). Acta Cryst. A52, 606–613. CrossRef CAS Web of Science IUCr Journals
© International Union of Crystallography. Prior permission is not required to reproduce short quotations, tables and figures from this article, provided the original authors and source are cited. For more information, click here.