Statistical
Interpolation of Sampled Images
ZHANG Yan, LI Xiaolong
1. Introduction
In practice, we use CCD, scanner or other devices to get the digital images, which have a finite resolution by sampling the images that have a high resolution. As a result, the sampled data are often aliased and may require interpolation. Many methods have been discussed for signal interpolation, such as bilinear interpolators, sinc and spline interpolators. However, a common property of all these existing methods is that the sampling process is not included in the design of the interpolators. In this paper, we will give an optimal Bayesian interpolator according to the idea presented in [1], which incorporates the sampling process and the statistics of the data and the noise.
The remainder of the paper is structured as follows. Section 2 discusses the sampling process, and the optimal interpolator is presented in Section 3. The simulation results are shown in Section 4. Then we draw the conclusion in Section 5.
All
real-world images are continuous. However, when they are acquired by image-capturing
devices, they must be represented in a digital form. Let us assume that the set
of samples representing the high-resolution data is
, where
When the same
scene is sampled at a lower resolution, the data is denoted by
, where
Here
where R
is the ratio between the sampling rates of
and
. Most image-capturing instruments adopt area sampling, where
a CCD cell integrates over an area to produce a sample. In the discussion that
follows, R is an integer greater than 1.
The
relationship between the high- and low-resolution data can be expressed as![]()
Where f, g and n are lexicographically ordered vectors of the unknown high-resolution data, the measured low-resolution samples, and additive noise values, respectively. All of these vectors are assumed to be stationary random processes. The matrix Q describes the sampling process. When the data-acquisition process of the image-capturing instrument is taken into account, Equation (1) ensures that the high-resolution samples are consistent with the low-resolution data gathered. We see that since R can be made arbitrarily high, the reconstructed data can approach continuous data.
There
are two possible sampling processes: point sampling and area sampling. We
assume
then consider a
4-neighbor structure:
![]()
Then we will use one point
to replace the
four points. For the point sampling, the sampling scheme is
; and for the area sampling, the sampling scheme is
. So, we can get
the sampling matrix easily. For example, assuming that the high-resolution data
are of size
(i.e. the length
of vector f is 64), then area sampling can be represented by

where each block
is

Similarly, we can get the point sampling matrix
for a situation
with the same dimensions. In point sampling process, we only use the
information of one point, which leads to lose much information; and the
simulation also show that the interpolation results are always worse than the
case of area sampling. So we confine our investigation to area sampling.
3.
Statistical Optimal Interpolation
3.1
Introduction to Optimal Interpolator
Following (1),
we consider the linear estimate
of
as
![]()
Where
is the optimal
interpolator that we want to find.
To find the optimal estimation, we view it as a Bayes
estimation problem associated with (1). Assuming that both
and
are multivariate
Gaussian processes, their respective probabilities are
![]()
![]()
Where
, which is the autocovariance of the high-resolution data;
and the autocovariance of noise on the measured low-resolution data is
. The data and the noise are assumed to be uncorrelated, i.e.
. Then we will find
such that the a
posteriori probability is maximized according to the Bayes rule:
![]()
Then the maximization problem is:
Whose solution leads to the optimal interpolator[2]
Here,
and
should be
invertible to ensure our discussion and the existence of the solution shown in (3). In fact, we know these two matrices are symmetric, so they
are positive definite. Hence the function in (2) that we
want to maximize is a strict concave function, so its solution is unique. In
addition, the bracketed factor in (3) is positive definite
and therefore invertible.
The sampling process and the statistics
of the data and the noise are incorporated in the optimal interpolator, which
leads to two main advantages of this optimal interpolator. Firstly,
describes the
statistical property of data. So, since different frequency components in the
high-resolution data lead to different frequency components in the
autocovariance matrix
, a better interpolation result can be obtained. Secondly,
the noise autocovariance matrix
works as a
balance between the fidelity of the interpolated image to the measured data and
the amount of noise present on the data.
3.2 Some
Assumptions of Autocovariance for High-Resolution Data
The above optimal interpolator relies on the knowledge of
the autocovariance
of the
high-resolution samples. This information may not be available with real world
data, but we can give some assumptions that are more realistic than others.
Every column and every row of an image can be considered as
a 1-D signal, whose autocovariance can be easily calculated (we note it as
) in the classical ways. Here, when we interpolate the 2-D image data, we take it as
a lexicographically ordered vector. So we can use
to get its
autocovariance:
,![]()
Then we give some classical assumptions
of
and show their
curves in Figure 1.
1) Gaussian model.
![]()
Which is consistent
with a smooth and continuous image. But some researchers suggest that the
Gaussian model may suffer from numerical instabilities and artifacts in
interpolation, and we think that it descends very quickly, so it can¡¯t estimate
the relationship of two points which are a little far from one to the other. So
we add a term to
to solve this
problem:
![]()
Here we take b=0.2 according to the experiments results.
2) Markov model.
![]()
Where c is the correlation between adjacent pixels.

Figure
1. Curves of autocovariance function
. Left: Gaussian model for a=2; Right: Markov
model for c=0.5.
3.3 Autocovariance
of Noise
The autocovariance of the noise is dependent on the correlation between the noise samples. There are two noise processes associated with image detection with a photoelectronic system: random thermal noise sources in the electronic circuits, and random fluctuations of the number of photons on the detector surface. We assume that the overall noise is uncorrelated.
The autocovariance of the noise is also a function of the
variance of the quantization process used in digital processing. In the case of
uniform quantization, it can be assumed that the CCD sampler has 8 bits, and
the images are captured on a gray scale with 256 levels. This leads to a
quantization step size
of 1, with the
variance of the noise being
. Since the noise is uncorrelated, the autocovariance is
![]()
Where
is the identity
matrix. In practices, the quantization procedure is often nonuniform. In this
case, the variance is larger at locations with more photons, and smaller at
locations with fewer photons. This means that the optimal interpolator can no
longer be applied globally to the data; instead, a different interpolator is
required for each point of interpolation. This increases the computational
complexity substantially. So in our paper, we confine the investigation to the
case of uniform quantization.
We use several high-resolution images (shown in Figure 2-(a), Figure 3-(a) and Figure 4-(a)) to test the optimal interpolator. Assume
, then we get the corresponding low-resolution images. Each
high-resolution image is divided into some
subimages, and each low-resolution image is divided into some
subimages. Then
each low-resolution subimage is interpolated by the optimal interpolator under
different autocovariance assumptions. Here, we take the Gaussian model for a=2.646
and the Markov model for c=0.98. And we use the normalized MSE defined
by (4) to evaluate the interpolation results. The results
are shown in Figure 2, 3 and 4. And the MSE of different models are shown in Table 1.
Table 1. MSE when using the classical interpolators (Bilinear and Cubic) and our optimal interpolators (Gaussian model and Markov model).
|
|
Bilinear |
Cubic |
Gaussian |
Markov |
|
Lena |
1.816 |
1.738 |
1.718 |
1.684 |
|
Goldhill |
8.405 |
7.633 |
7.205 |
6.412 |
|
Room |
8.055 |
6.357 |
3.486 |
3.917 |
From these images, we see that, in general, the MSE of the optimal interpolated results is smaller than the ones of the cubic and bilinear interpolated result. Particularly, in our optimal interpolators, the Markov model performs better than Gaussian model. In addition, the bilinear and cubic interpolation results are a little blurry while our interpolation results are more net. When interpolating the textured or detailed image areas, the optimal interpolator shows an advantage. For example, for image Room (Figure 4(d)), the stripes on the towel are reconstructed well using our interpolator; and the results are blurry when using classical interpolators. So we think the idea of using statistical model to design an optimal interpolator is a meaningful method.
Since we interpolate the image block by block, the block-effect may appear. But this effect can be avoided by choosing the proper parameters and by increasing the size of block.
Now we add the term
and try to
denoise. Figure 5 shows the result of interpolation when
taking
. The MSE is a little smaller. In Figure 6,
we add random noises to image Lena (Figure 6-(a)) so that
it can¡¯t be recognized. Then we use Markov model and take
to realize the
interpolation. We see that cubic interpolator can¡¯t reduce noise, but our
optimal interpolator can reduce some noises and make the noised image a little
recognizable, and we can expect a much better result when we increase the size
of block.

(a) Lena:
original
(b)
Cubic Interpolation (MSE=1.738
).

(c) Gaussian model (MSE=1.718
). (d) Markov model (MSE=1.684
).
Figure 2. Different interpolation of Lena image. (b) use classical cubic interpolator; (c) and (d) use our optimal interpolators. The MSE is given for different interpolators.

(a) Goldhill: original
(b) Cubic Interpolation
(MSE=7.633
).

(c) Gaussian model (MSE=7.205
).
(d) Markov model (MSE=6.412
).
Figure 3. Different interpolation of Goldhill image. (b) use classical cubic interpolator; (c) and (d) use our optimal interpolators. The MSE is given for different interpolators.

(b)
Cubic Interpolation (MSE=6.357
).

(c)
Bilinear Interpolation (MSE=8.055
).
(d) Optimal interpolation using Gaussian model (MSE=3.486
).

(e) Optimal interpolation using Markov model (MSE=3.917
).
Figure 4. Different interpolation of Room image. (b) and (c) use classical cubic and bilinear interpolators; (d) and (e) use our optimal interpolators. The MSE is given for different interpolators.

Figure 5. Optimal interpolation of Lena image
using Markov model and taking
. The MSE is 1.683
.


(b) Cubic interpolation (c) Optimal interpolation
Figure 6. Optimal interpolation of noised Lena
image using Markov model and taking
.
We realize an optimal interpolation method according to the ideas given in [1]. Compared to classical cubic and bilinear interpolators, our optimal interpolator performs better, but the computational complexity is increased since it needs to calculate the inverse of a large matrix. In addition, this interpolator incorporates the sampling process and the statistics of noise, which ensure that our interpolator can reduce some noises when interpolating. This interpolator depends on the choice of the parameter of the autocovariance and as a result our choices in this paper may not adapt to all images. Thus, how to select a parameter adaptively and quickly is a problem worthy to investigate.
Remarks:
1) In paper [1], when calculating the MSE, the authors removed the block edges to eliminate the influence of edge effects. But we think it¡¯s not a reasonable evaluation criterion. So here we use all data to calculate the MSE, and when using our chosen parameters, there is few block effects.
2) The authors of paper [1] demonstrated that this interpolator can partially remove the aliasing. But we think this is not the case in general. We try it to several images and find that the aliasing still exist.
References:
[1] W.-Y.V.Leung, P.J.Bones and R.G.Lane, ¡°Statistical interpolation of sampled images.¡± in Optical Engineering, 40(4), pp547-553, Avril 2001.
[2] H.C.Andrews and B.R.Hunt, ¡°Digital image restoration.¡± Prentice-Hall Signal Processing Series, Prentice-Hall, Englewood Cliffs, NJ, 1977.
Source Code
samplematrix.m : generate the sampling matrix ![]()
sample.m : generate the low-resolution data
kernel.m :
generate the autocovariance matrix of data and calculate
defined in (3)
main.m : main interpolation algorithm