IsmailPeters

Back to Psych 221 Projects 2013

Background

Camera Basics

A digital camera captures color images by recording data in different channels of the visible spectrum (red, green, blue). Most cameras accomplish this by using a single type of CCD or CMOS sensor at each pixel. This sensor measures the intensity of light that is focused on it, but cannot distinguish between colors. Therefore, a color filter is placed over each sensor, and the sensor only records the intensity of a specific color at that pixel. The arrangement of these color filters over an image is known as a Color Filter Array. Figure 1 shows a very common CFA known as the Bayer Array.

Figure 1: Bayer Array

In the Bayer Array, green color filters are arranged in a checkerboard pattern, whereas the red and blue color filters are arranged by alternating rows. The original reasoning for using this distribution of colors was to mimic the physiology of the human eye. During daylight vision, the luminance perception of the retina uses both L and M cones, which are more sensitive to green light. The red and blue filters control the sensitivity of the eye to chrominance.

CFA Interpolation

Since each color channel is only tallied at specific coordinates, the remaining pixels in that channel must be estimated in some way. There are many interpolation schemes - a few simple examples are the bilinear,bicubic, and smooth hue transition interpolations shown below. Other interpolation schemes, such as median filter, gradient-based, and adaptive color plane, become increasingly complex. A wide survey of interpolation methods is available, along with an analysis of the advantages and disadvantages of each.

Bilinear/Bicubic

A linear combination of the nearest N neighbors in any direction for N = 1 (bilinear) or N = 3 (bicubic). For a given matrix of known values R(x,y), the interpolated set of values r(x,y) can be computed:

$r(x,y)=\sum _{u,v=-N}^{N}h_{r}(u,v)R(x-u,y-v)$

where h_r is the 2-D kernel containing interpolation weights based on relative neighbors of a given pixel.

$h_{1}={\frac {1}{4}}{\begin{bmatrix}0&1&0\\1&4&1\\0&1&0\end{bmatrix}},h_{3}={\frac {1}{256}}{\begin{bmatrix}0&0&0&1&0&0&0\\0&0&-9&0&-9&0&0\\0&-9&0&81&0&-9&0\\1&0&81&256&81&0&1\\0&-9&0&81&0&-9&0\\0&0&-9&0&-9&0&0\\0&0&0&1&0&0&0\end{bmatrix}}$

Smooth Hue Transition

Operates on the assumption that the hue of a natural image (chrominance/luminance, or red/green and blue/green) changes smoothly in local regions of that image. In order to calculate this, the green channel is first interpolated using bilinear interpolation. A missing red point r(x,y) is calculated from known red values R(x,y) and green values G(x,y).

$r(x,y)=G(x,y)*{\frac {1}{2}}\left({\frac {R(x-1,y)}{G(x-1,y)}}+{\frac {R(x+1,y)}{G(x+1,y)}}\right)$

The sum shown above is valid for when the interpolation pixel has two adjacent pixels in the same row. The equation is altered slightly for two adjacent pixels in the same column, or when there are four adjacent pixels in the corners. The missing blue pixels are calculated in the same way.

Detecting Forgeries

Regardless of the type of interpolation used in a given image, the estimated pixels should exhibit a strong correlation, or dependence, on their surrounding pixels. If one can categorize each pixel as either correlated or independent with respect to other pixels, one should see a periodic pattern that mimics the CFA used to construct the image. Even if an altered image does not have any visual cues that point to its forgery, inspecting the correlation of pixels to one another can potentially expose which parts of an image were altered. It is important to note that, even if an image has been altered, the mentioned correlations may still be kept intact. In this way, detecting forgeries through CFA interpolation is one tool out of many used to detect forgeries.

Methods

Finding Correlations

The first step in detecting forgeries is categorizing the pixels in an image into two groups: those that are correlated with neighboring pixels and those that independent. Since the camera, CFA, an interpolation method are unknown for a given image, a robust algorithm must be used to iteratively deduce this categorization.

Expectation/Maximization

The EM algorithm has a wide variety of uses for working on a statistical set of data. In the context of this problem, each step of the algorithm has a specific purpose. In the expectation step, the probability of each pixel in a given image being correlated to its neighbors is calculated based on the latest guess of interpolation coefficients. In the maximization step, new interpolation coefficients are calculated based on the probability map that reveals which pixels are indeed interpolated. Iterating through this process converges to a stable estimate of both the probability map and interpolation coefficients. If a pixel is correlated to neighboring pixels, it takes the form:

$f(x,y)=\sum _{u,v=-N}^{N}\alpha _{u,v}f(x+u,y+v)+n(x,y)$

where $\alpha _{u,v}$ are the interpolation coefficients and n(x,y) is the residual error of the linear model and ideally a Gaussian distribution with zero mean. The probability that a sample f(x,y) exists and that it is correlated to its neighbors is given by the equation:

$P(x,y)={\frac {1}{\sigma {\sqrt {2\pi }}}}exp\left[-{\frac {1}{2\sigma ^{2}}}\left(f(x,y)-\sum _{u,v=-N}^{N}\alpha _{u,v}f(x+u,y+v)\right)^{2}\right]$

where $\sigma$ is the standard deviation of this Gaussian distribution and is updated in the M-step. The new estimate of the interpolation coefficients is found by minimizing the quadratic error function:

$E(\alpha )=\sum _{x,y}w(x,y)\left(f(x,y)-\sum _{u,v=-N}^{N}\alpha _{u,v}f(x+u,y+v)\right)^{2}$

By setting the derivative of this equation equal to zero and rearranging terms, the following equation exists for each relative index (s,t) of the set interpolation coefficients $\alpha (u,v)$:

$\sum _{u,v=-N}^{N}\alpha _{u,v}\left(\sum _{x,y}w(x,y)f(x+s,y+t)f(x+u,y+v)\right)=\sum _{x,y}w(x,y)f(x+s,y+t)f(x,y)$

where w(x,y) = P(x,y). By performing this for each coefficient, a system of M equations and M unknowns, M being the number of interpolation coefficients, is set up. In this way, a new series of coefficients can be computed.
In order to make the computation optimized for MATLAB, the summations were converted into matrix multiplications. For example, for a given interpolation coefficient $\alpha _{s,t}$ , the w(x,y) , f(x+s,y+t) , and f(x,y) matrices can be written as row vectors of size 1 x N where N is the number of pixels under analysis. Thus the right hand side of the equation can be represented with matrix multiplations:

$\sum _{x,y}w(x,y)f(x+s,y+t)f(x,y)={\begin{bmatrix}w(x,y)\end{bmatrix}}.*{\begin{bmatrix}f(x+s,y+t)\end{bmatrix}}*{\begin{bmatrix}f(x,y)\end{bmatrix}}^{T}$

Similarly, the left hand side of the equation can be re-written. If all interpolation coefficients are taken into account, the entire system of equations can be represented:

${\begin{bmatrix}w(x,y)\\.\\.\\w(x,y)\end{bmatrix}}.*{\begin{bmatrix}f(x+1,y+0)\\f(x+0,y+1)\\.\\f(x+s,y+t)\end{bmatrix}}*{\begin{bmatrix}f(x+1,y+0)\\f(x+0,y+1)\\.\\f(x+u,y+v)\end{bmatrix}}^{T}*{\begin{bmatrix}\alpha _{1,0}\\\alpha _{0,1}\\.\\\alpha _{u,v}\end{bmatrix}}={\begin{bmatrix}w(x,y)\\.\\.\\w(x,y)\end{bmatrix}}.*{\begin{bmatrix}f(x+1,y+0)\\f(x+0,y+1)\\.\\f(x+s,y+t)\end{bmatrix}}*{\begin{bmatrix}f(x,y)\\f(x,y)\\.\\f(x,y)\end{bmatrix}}^{T}$

For clarity, the resulting sizes of these matrices is contained below.

${\begin{bmatrix}MxM\end{bmatrix}}*{\begin{bmatrix}\alpha _{1}\\\alpha _{2}\\.\\\alpha _{M}\end{bmatrix}}={\begin{bmatrix}Mx1\end{bmatrix}}$

In this way, a new set of coefficients can then be easily computed. This rearrangement of the data provided much faster computations, and it was very easy to add more coefficients to the system. The EM algorithm was allowed to run until the coefficients generated in successive iterations differed by an amount less than some arbitrarily small threshold.

2-D FFT

Once a probability map was generated from the EM algorithm, the next step was to determine if the periodic pattern of the CFA was present in the map. To do this, a two-dimensional FFT was performed on the probability map. The default 2-D FFT magnitude plot in MATLAB contains all of the low frequency information in the upper left corner, and the high frequency components along the left and bottom edges. The plot also scales to the largest value. In order to make the desired periodicity easier to observe, the FFT plot was up-sampled, high-pass filtered, and shifted so that the low frequency components are at the center of the image. Figure 2 shows the original FFT plot for a given probability map. Figure 3 is the same data, but it has been upsampled and filtered in order to bring forward the periodic effects of the CFA interpolation.

Figure 2: Raw FFT

Figure 3: Improved FFT

Similarity Measure

In order to build a classifier, a reliable quantity must be measured that is different by some threshold for real and fake images. The method chosen utilizes the similarities between the FFT of the probability map of an image and the FFT of an ideal synthetic grid. To generate a synthetic grid for a given color channel, any pixel that corresponds to that channel's color is marked with a 1 at that location; otherwise, all other locations contain a zero. In order to avoid the event that the probability map and synthetic grid are out of phase, which would yield a low similarity when compared directly, the magnitudes of the FFT results are compared. If $P(\omega _{x},\omega _{y})$ is the FFT of the probability map, and $S(\omega _{x},\omega _{y})$ is the FFT of the synthetic grid, then the similarity measurement M is:

$M=\sum _{\omega _{x},\omega _{y}}|P(\omega _{x},\omega _{y})|\cdot |S(\omega _{x},\omega _{y})|$

This measurement is essentially a dot product of the FFT data for a syntehtic image and the image under consideration. For windows in which the probability map is close to that expected, this similarity measure should remain constant. If the probability map distribution varies significantly, the similarity measure will deviate from the average value. This can be true for any channel. Figure 4 shows the FFT of the generated synthetic grid.

Figure 4: Synthetic Grid FFT

Results

The test set used to examine the performance of the work described above consisted of 70 image pairs. In each image pair, one image was contained an unaltered, CFA interpolated, 'real' image. The other image in the pair was visually identical, but had undergone some alterations. A window size of 128 pixels x 128 pixels was used to examine these images.

EM Algorithm Performance

For a majority of windows in all of the images, the EM algorithm converged to a solution within 100 iterations. The most troubling cases were windows of uniform color - it was difficult to mathematically converge to a solution similar to other windows given the uniformity of the window. The other cases that did not easily converge were those that contained altered parts of an image. Figure 5 shows a real image with the test window marked; the probability map generated by the EM algorithm, and the FFT plot of the probability map.Figure 6 shows the same information for the fake image in the same image pair.

Figure 5: Real Image, Probability Map & FFT

Figure 6: Fake Image, Probability Map & FFT

Comparing Real and Fake Image Pairs

After seeing that both the real and fake images show very periodic results for portions of the image pairs that are the same, the similarity measure was examined for the image pairs. Figures 7 through 12 show the similarity measurement for a 128 x 128 pixel window that was swept in 1/2 window step sizes for several images.

Figure 7: Similarity of Real Image Probability Map to Ideal

Figure 8: Similarity of Fake Image Probability Map to Ideal

Figure 9: Similarity of Real Image Probability Map to Ideal

Figure 10: Similarity of Fake Image Probability Map to Ideal

Figure 11: Similarity of Real Image Probability Map to Ideal

Figure 12: Similarity of Fake Image Probability Map to Ideal

Setting a Threshold

From the previous set of figures, it is clear that some image pairs produce significant differences that are separated by a clear threshold. However, for other images, the differences are not so clear cut. For the entire test set, the distribution of the similarity calculation was so large, it was impossible to determine a single threshold to distinguish between real and fake images. This can be seen within the following figures. One of the techniques used to increase the accuracy of the classifier was to exclude irregularities unless they occur in all three color channels at the same location. If one color channel varied significantly, but the other two exhibited normal behavior, an image was not considered to be 'tampered'. Figures 13 and 14 show the problem with finding a threshold that works for all images; the real image gives mixed results due to the almost constant background color.

Figure 13: Fake Image with Alleged Tamper Window

Figure 14: Real Image with Falsely Identified Tamper Window

JPEG Compression

The classifier was tested for varying levels of JPEG compression to examine the robustness. For an image that was JPEG compressed for qualities ranging from 10 - 100, the altered portion of the image was visible until the compression level reached 50% as seen below. For all of the JPEG compressed images, the algorithm to detect tampering through CFA interpolation still correctly identified doctored images for quality factors as low as 50%.

Figure 15: Series of images at different compression levels

Figure 16: Corresponding similarity plot

Final Test Images

Of the final 11 sets of test images that were provided, the classifier was able to correctly categorize 5 properly(fake classified as fake, real classified as real). The algorithm is still highly sensitive to re-sizing of the entire image and JPEG compression after the image has been modified. Figures X,Y demonstrate the issues seen with jpg compression and why it is difficult to determine a single threshold value for both images. Similarly, figures 19 and 20 demonstrate the issues seen by re sizing the tampered image.

Figure 17: Similarity plot for tampered png image

Figure 18: Similarity plot for tampered jpg image

Figure 17: Similarity plot for tampered png image

Figure 18: Similarity plot for re-sized png image

Conclusions

The algorithms provided could still use further improvements as they stand. First, the EM algorithm has trouble working with image windows of constant color. This low spatial frequency content make it hard for the algorithm to separate the pixels correlated with their neighbors and the ones that are independent. If this detail can be teased out of the algorithm, the similarity plots would be much more stable for many more of the images. Another improvement that can be made is figuring out how to increase the accuracy of the classifier. As of right now, the classifier is very sensitive to re-sizing or JPEG compression of the image after it has been modified. Possible methods to improve the hit rate would be decreasing the image window and step size. One possible reason that the accuracy is reduced is that the sweep window is too large and not sensitive to tampering.
As it stands, the classifier can successfully identify a fake image as fake and real image as real for 41 test pairs of the original 70 test pairs given.

References

1. Popescu, A.C.; Farid, H.,Exposing digital forgeries in color filter array interpolated images, Signal Processing, IEEE Transactions on , vol.53, no.10, pp.3948,3959, Oct. 2005

2. R. Ramanath, W. E. Snyder, G. L. Bilbro, and W. A. Sander III, “Demosaicking methods for Bayer color arrays,” Journal of Electronic Imaging, vol. 11, no. 3, pp. 306–315, July 2002.

3. A. Dempster, N. Laird, and D. Rubin,“Maximum likelihood from incomplete data via the EM algorithm,” Journal of the Royal Statistical Society, vol. 99, no. 1, pp. 1–38, 1977.

Appendix I - MATLAB Code

Code

File:IsmailPeters CameraForensics.zip Source Code

Appendix II - Work partition

Partners through and through. 50/50.
Although Chasen developed all of the sweet UI tests.