Reconstructing Images From Their Most Singular Fractal Manifold

Real-world images are complex objects, difficult to describe but at the same time possessing a high degree of redundancy. A very recent study on the statistical properties of natural images reveals that natural images can be viewed through different partitions which are essentially fractal in nature. One particular fractal component, related to the most singular (sharpest)transitions in the image, seems to be highly informative about the whole scene. In this paper we will show how to decompose the image into their fractal components.We will see that the most singular component is related to (but not coincident with) the edges of the objects present in the scenes. We will propose a new, simple method to reconstruct the image with information contained in that most informative component.We will see that the quality of the reconstruction is strongly dependent on the capability to extract the relevant edges in the determination of the most singular set. We will discuss the results from the perspective of coding, proposing this method as a starting point for future developments.


Reconstructing Images From Their Most Singular Fractal Manifold Antonio Turiel and Angela del Pozo
Abstract-Real-world images are complex objects, difficult to describe but at the same time possessing a high degree of redundancy.A very recent study [1] on the statistical properties of natural images reveals that natural images can be viewed through different partitions which are essentially fractal in nature.One particular fractal component, related to the most singular (sharpest) transitions in the image, seems to be highly informative about the whole scene.In this paper we will show how to decompose the image into their fractal components.We will see that the most singular component is related to (but not coincident with) the edges of the objects present in the scenes.We will propose a new, simple method to reconstruct the image with information contained in that most informative component.We will see that the quality of the reconstruction is strongly dependent on the capability to extract the relevant edges in the determination of the most singular set.We will discuss the results from the perspective of coding, proposing this method as a starting point for future developments.

I. INTRODUCTION
E DGE detection is a common feature of the mammals' vi- sual neural system [2], [3].It has been proposed that edge detectors could be used to provide efficient coding algorithms [4] and in fact maximization of the information transfer lead to orientational edge-detecting filters [5].However, providing a reasonable, nonconventional definition of "edge" is more controversial [6].
A different strategy to produce efficient coding can be that of the statistical analysis of images [7], [8].This kind of analysis implies to identify the origin of redundancies (like that of the power spectrum [9]) to devise redundancy-reducing codes (as in [10], [11] for the case of the power spectrum).
In the last years a new statistical study has arisen, that of multifractals in natural images [1], [12].The multifractal scheme provides a richer framework than that of the simple characterization of the power spectrum.This scheme makes possible to split any image into a collection of fractal sets, from which one of them is supposed to be the most informative [1].Not surprisingly, that set is usually edge-like [1], [13].
In this paper we will make use of the multifractal scheme to obtain the most informative component.We will propose a method to reconstruct all the image just using information contained in this set.We will design this method by requiring it to have several reasonable features.We will study its theoretical properties and discuss its experimental performance over real data.
The paper is organized as follows: In the following Section some notations and the methods are defined.Section III is devoted to review the multifractal scheme.In Section IV the reconstruction procedure is issued and in Section V the quality of the method over real natural images is tested and discussed.The statistics about the computer implementation of the new techniques issued are presented in Section VI.In Section VII we discuss the results and some of the possible improvements.Finally, the conclusions are presented in Section VIII

II. NOTATION AND METHODS
We will denote the recorded field of luminance intensities as .We will work with an additive normalization of it, the global contrast field, which is defined as (1) where is the average luminosity over the image.So, the average of over any image vanishes.As monochrome test images we will make use of several pictures taken from Hans van Hateren's web database (see [14] for technical details).In all the cases we will work on several 512 512 patches from those scenes.We will also use the 512 512 resolution greylevel version of Lena's picture to illustrate several examples.
The determination of the multifractal structure of each image was made by means of wavelet analysis [15], [16] using wavelets from the family [1] for and averaging the resulting coefficients.The distribution of exponents was computed for each image recording the relative frequencies.The most singular exponent was computed as the mean of the 1% and 5% quantils.The dispersion around this value is conventionally fixed in a convenient value, depending on the image.

III. FRACTAL DECOMPOSITION OF IMAGES
There are several equivalent ways to show that natural images possess multifractal structure.One possibility consists of constructing a positive measure , which assigns a positive value to any set .The measure can be defined by its density as: (2) We will define a measure which takes into account the sharp transitions found inside each area.Following [1], we will define the measure density as (3) The measure gives an idea of the local distribution of the gradient of and its inhomogeneities across the image.It is possible to use it to characterize the behavior of any particular point by observing the evolution of the measure of balls centered on .By we will denote the ball of radius centered on .The measure is said to define a multifractal if at each point in the image, the measure can be characterized by a local singularity exponent in the way (4) where means a term which is negligible in comparison with for small values of and is the dimension of the space.For a multifractal measure, (4) determines uniquely the coefficient and the exponent : they can be obtained by linear regression of versus .The coefficient depends on some arbitrary choices, like the particular metrics used to define the balls and the scale unit for and it provides no information about the changes in scale.On the contrary, the exponent is independent of the metrics and scaling unit and gives all the information about the evolution under changes in scale (changes in ; see [1] for a full discussion of the interpretation of ).We will see that there exist better choices than to complement the information provided by (for instance, the value of over a particular subset).
Natural images have been shown to exhibit multifractal behavior with the measure defined in (3) (see [1], [13]).However, in practice, a direct application of log-log linear regression on (4) yields rather coarse discrimination of the exponents, specially due to discretization of the values in the radius for the balls (see [17] for a discussion about how to deal with discrete circles).The standard technique used to circumvent this difficulty involves the use of wavelet projections as singularity analyzers.Let be an appropriate function and let us define the wavelet projection of over at the point and scale as (5) It can be proven [15], [16] that for multifractal measures (i.e., those verifying (4)) the wavelet projections also scale as power laws in , namely (6) where is the same as in (4) but depends on the choice of .The exponents are obtained as the slope of the linear regression of versus .This method has very good performance in practice [1], [13].
The existence of a multifractal measure implies a strong hierarchical organization in images.This "multiple-fractal" character shows up when the image is split in the different singular components , being these formed by the points sharing the same singularity exponent (7) The components are observed to have nontrivial fractal dimensions [1], which can be predicted from statistical properties of the images [18].For that reason, those sets are usually called the fractal components of the image (see Fig. 1).Over real images, a direct numerical estimation of for each fractal component is possible once the multifractal decomposition is performed.However, some popular techniques as that of box-counting dimension (see [19] and references therein) lead to overestimated values of the dimension because fractal components are usually topologically dense.For that reason, in many contexts [20] it is preferred to analyze the statistical structure with wavelet-based techniques to obtained the correct .Fortunately, for wide collections of natural images the observed multifractals belong to the class of the log-Poisson multifractals [1], [12], [21].The function (usually called singularity spectrum) of log-Poisson multifractals can be described in terms of only two free parameters and so no complicated techniques are necessary in its determination.The general expression for a log-Poisson singularity spectrum is given by [1] (8) where . The two required free parameters are the minimum possible exponent and the dimension of the fractal component associated to it, [1], [12].That fractal component, denoted , is called the most singular manifold (MSM).
The MSM plays a fundamental role in the context of multifractals in natural images.It is usually observed [1], [12], [13] to be of dimension 1.0 but with noninteger singularity exponent.Visual inspection of this set reveals a structure which resembles the "edges" or contours present in the scene [1], [13], which would take account of the dimensionality 1.0 of the set (Lena's image has also ; see Fig. 1).For that reason, we will always assume that , so to define completely the multifractal we will just need to calculate (which is estimated as explained in Section II).
As log-Poisson multifractals can be described in terms of the parameters defining the MSM, it was proposed in [1] that the MSM could contain the majority of or all the information conveyed by the image.This would imply, for instance, that the second and third manifolds represented in bottom of Fig. 1 could in fact be calculated from the first (top right).We will see in the next section that it is possible to propose a propagator for reconstructing the images from information which is contained in the MSM.

IV. DETERMINISTIC RECONSTRUCTION
We will propose a simple, deterministic propagator for the multifractal measure from its restriction to the MSM.So, we should consider the measure density and reconstruct it  from its values over .The choice of the propagator is made under the requirement of verifying the following five conditions: 1) it is deterministic; 2) it is linear; 3) it is translationally invariant; 4) it is isotropic; 5) it leads to the observed power spectrum.Under this set of assumptions, which we will next explain, there exists only one possible propagator.This theoretical propagator (that we will also call "reconstructor") will or will not reconstruct the image from the data, which should be experimentally checked.Now, we will construct the propagator following the properties we have required.
The deterministic character of the reconstruction allows to consider the propagator not as a random variable at each point, but as an actual function of its arguments.The following expression holds: (9) that is, the value of at any point is a function of the values of over the MSM.The linearity of the propagator forces us to work with the gradient of the contrast, , instead of its modulus (the modulus is a nonlinear operator).Anyway, the gradient itself contains more information than its modulus, so if the reconstruction is possible from the modulus, it is also possible from the gradient.Hence it follows that there exists an integral representation for the function in the way (10) where means line integration along the MSM.We are representing the linear operator by means of its density, denoted .That density is a 2 2 matrix; we will represent the matrix element as where .If we denote by , the components of , the vector can be represented also by its coordinates where (11) The vector represents the vectorial density of the gradient, because when integrated over the MSM it turns out the value .The translational invariance is a usual requirement in image statistics, meaning that there is not a preferred place in which objects could be expected to be found at natural scenes.Its experimental extent is limited due to the finite size of images.In terms of the integral representation, (10), it implies that the operator density is in fact a function of , that is (12) In fact this equation can be simplified to a scalar equivalent: the left hand side is a gradient (perhaps in a distributional sense, as the contrast is discontinuous), so its curl vanishes: .This property is translated to the reconstructor , in the way (13) which implies (14) that is, is the gradient of a vector .So ( 12) can be simply expressed as (15) The above equation can be rewritten in a very useful form defining the field as (16) where stands for the density of the proper Hausdorff measure restricted to the set (a delta function over the lines defining ).In this way, (15) becomes a convolution, because now the integration is performed over all the space and no longer over the lines of the MSM (the restriction is indeed still present, but now it is introduced by the field ), that is (17) where stands for the convolution operator.So, the reconstruction formula is elegantly expressed in the Fourier space as (18) which is an integral equation equivalent to (15); let us notice that " " means the scalar product of the complex vectors.Recall that now the boundary conditions are contained in the vector field , which depends on the particular image to be reconstructed.The reconstructor is defined by the complex vector field .It is thus natural to require the reconstructor to be isotropic, as the particularities of the image are already contained in and we think that is an universal propagator.This implies that (19) where stands for the modulus of the propagator.Because of the isotropy, can only be a function of the modulus of the frequency vector.
To define completely the propagator, we recall a well established property of natural images, namely the scaling of their power spectrum (see [9]), which is (20) where is a nonuniversal, small exponent which depends on the particular image ensemble considered (see for instance [22]).The simplest possible is then given by , is where by we denote the imaginary unit, .The modulus of the Fourier transform of the contrast is then given by , where has a weak dependence on and varies from one image to another.Hence, according to the definition of the power spectrum (22) so the term gives rise to the factor while the factor introduces the particular anisotropies of the image and would give rise to the weak dependence in (20).So finally (21) defines our propagator and the reconstruction formula, (18), reads (23) Practical application of (23) on discretized images is quite simple.The procedure is as follows.
1) The singularity exponents at each point on the image are computed as in [1].
2) From the distribution of singularity exponents the value of the most singular exponent is calculated as the average of the 1% and 5% quantils.The dispersion around this value is conventionally fixed ( 0.2 is usually a good choice).
3) We define the density function as 1 if the exponent associated to the point can be considered as the most singular one (i.e., ) and 0 otherwise (i.e., or ).4) We compute the gradient and we obtain the essential vector field as (that is, it equals the gradient over the MSM and it vanishes outside the MSM) 5) The bidimensional vector field is Fourier transformed to obtain the complex bidimensional vector field 6) The scalar product with the frequency vector is computed: 7) The complex number so obtained is multiplied by the imaginary unit and divided by to obtain .An anti-Fourier transform provides the reconstructed .The reconstruction formula has a very interesting property: no matter the image considered, (23) allows reconstructing the correct provided that the set (which defines ) is large enough: if is taken as the whole image, and (23) turns out to be a trivial identity.The question is if natural images allows reconstruction considering a rather sparse set : the MSM.

V. EXPERIMENTAL PERFORMANCE
The experimental application of the reconstruction formula is very simple.First, the singularity exponents are computed for each point in the image analyzed.Then, the least exponent is found and the associated MSM is isolated.To compute the field , we compute the gradient ; according to (16) equals over the MSM and zero outside.The application of the reconstruction formula, (23), in the Fourier space is then straightforward.
In Fig. 2 we show the field which was obtained from a coarse version of the MSM.It is rather unclear which should be the appropriate quantization noise in and we have taken the one which allows a reasonable reconstruction out of the MSM; the central value , however, is well determined: it is fixed by the 1%-5% quantils average described in Section II.It could be argued that in fact there is no a sparse MSM capable of reconstructing via the reconstruction formula, but a collection of sets of increasing size obviously reconstruct better and better.However the problem seems to be the quality of the edge detection.In Fig. 3 the error image for Lena's picture is shown.The error image is defined as the difference between the original contrast and the one generated reconstructing from the MSM, .This error image seems rather fuzzy and the only recognizable features are, precisely, the edges that were not previously detected.
The performance of the reconstruction over Lena's image seems good, but the quality is probably lowered by the (unknown) filters applied on that particular image.This is reasonable because filtering damages the natural propagation of light, thus the reconstruction is likely to work much better on nonprocessed natural images.In this spirit we repeated the process with other, non-filtered images (see Fig. 4 for several examples).The general performance is good, although if one border is lost (at the time of edge detection) so it is all the structure associated to it, which seems quite reasonable.This problem could have also to do with fluctuations in the multifractal structure [21].
The quality of the reconstruction thus varies largely from one place to the other in the image: at those places in the neighborhood of a lost edge the differences between the reconstructed contrast and the original one, , are large; at the places contained inside of a well defined object, the differences are small.To give a reasonable measure of this inhomogeneus error we computed the ratio of the power spectrum of the error image to the power spectrum of the original contrast; the result is shown in Fig. 3.It seems that the error is inhomogenously distributed across the frequencies, being more important at lower and higher frequencies.This fact could correspond to raylight spreading effects caused by the loss of some edges (lower frequencies) and the loss of the edges itself (higher frequencies).

VI. COMPUTATION STATISTICS
We ran all our programs over a DEC WorkStation at 500 MHz working under UNIX and always for 512 512 images.All the programs are written in C but they are not optimized.The programs used to isolate the fractal components (the MSM in particular) took 27 s to produce the output files, with a peak of memory use of 13 Mb.The computation of the reconstruction data (the vector field ) took 10 s and 9 Mb of memory usage.The reconstruction out the vector field took an insignificant amount of time, under 3 s and a memory usage of 9 Mb.The MSMs were recorded as two-color images (e.g., black for the MSM, white for its complementary).A measure of the complexity of the MSM was given by its density, that is, the number of points belonging to the MSM with respect to the total number of points.We will call that density "edge density."The edge densities range from 10% for simple images to 40% for very complicated images, the values ranging from 25 to 30% being typical.According to the mentioned values of density, a naif coding of the vector is then almost not useful from the point of view of compression, as this vector involves two numbers (two coordinates) for location of the MSM (however, the vector varies smoothly along the edges and it is almost perpendicular to them; two facts we will use to devise better codings in next works) The PSNRs of the reconstructions ranged from 20 to 40 dB, the range 27-35 dB being typical.This value is not obviously related to those of the edge densities, being observed low edge density images with excellent PSNRs and high edge density images with poor PSNRs.We concluded that the quality of the reconstruction is apparently independent of the edge density.

VII. DISCUSSION
This description of the image out of the MSM seems to be, at least, a reasonable first order approximation.There are three possible interpretations for the observed deviations and so three possible ways to correct them.
1) The reconstruction formula is correct and the deviations are due to lack of accuracy in the computation of the MSM the computation of the MSM should be improved.
2) The reconstruction formula is correct, the MSM is correctly computed, but the MSM alone does not describe completely the image the other informative structures should be identified.
3) The reconstruction formula is not correct, but the MSM contains all the relevant information an improved reconstructor should be proposed.The first possibility implies that better techniques of edge detection should be devised to take all the potential out of this technique.The second possibility would probably mean a deviation from the Log-Poisson model, which is contradictory with the present experimental evidence (see [1], [12]), unless those deviations are very small but significant in the context of each image.The third possibility could be interesting to start a new program.

VIII. CONCLUSIONS
In this paper we have reviewed the multifractal formalism applied to real world natural images, which allows to split the image into a collection of fractal sets.These fractal components can be described using just one of them, namely the Most Singular Manifold (MSM).We have proposed a propagator to reconstruct the whole image using the gradient over the MSM.We have discussed the validity of this method and its quality over real data.Besides, we have proposed three different ways to improve the performance.
The relevance of the MSM as the most informative set in the image (already proposed in [1]) emphasizes the role of the contours and the objects in natural images, which could be connected to other coding algorithms [5] and with biological observations [2], [3].
The possible applications of the method presented here go from coding schemes for real world images to image processing and analysis.
Authorized licensed use limited to: Universitat de Barcelona.Downloaded on February 16, 2009 at 06:32 from IEEE Xplore.Restrictions apply.
Authorized licensed use limited to: Universitat de Barcelona.Downloaded on February 16, 2009 at 06:32 from IEEE Xplore.Restrictions apply.

Fig. 3 .
Fig.3.Left: Error image for Lena's picture.This image can also be generated from the complementary of the reconstructing manifold.The level of contrast of this image has been linearly increased to enhance the details.Right: Ratio power spectrum of the error image-power spectrum of the image, for Lena's image (continuous line) and its average over the 512 2 512 central patches of 1000 random van Hateren's images (dashed), in semi-log plot.The power spectra were radially averaged.The error concentrates in the higher and the lower frequencies.