KT: kaleidoscope transform
MRI: magnetic resonance imaging
DFT: discrete Fourier transform
CS: compressed sensing
ChaoS: chaotic sensing
CT: computed tomography
PET: positron emission tomography
RMSE: root-mean-squared-error
PSNR: peak signal-to-noise ratio
SSIM: structural similarity
VIF: visual information fidelity
SVD: singular value decomposition

Bespoke Fractal Sampling Patterns for Discrete Fourier Space via the Kaleidoscope Transform

Jacob M. White, , Stuart Crozier, and Shekhar S. Chandra J. M. White, S.S. Chandra and S. Crozier are with the University of Queensland, Brisbane, St Lucia QLD 4072, Australia (e-mail: [email protected])

Abstract

Sampling strategies are important for sparse imaging methodologies, especially those employing the discrete Fourier transform (DFT). \AclChaoS is one such methodology that employs deterministic, fractal sampling in conjunction with finite, iterative reconstruction schemes to form an image from limited samples. Using a sampling pattern constructed entirely from periodic lines in DFT space, chaotic sensing was found to outperform traditional compressed sensing for magnetic resonance imaging; however, only one such sampling pattern was presented and the reason for its fractal nature was not proven. Through the introduction of a novel image transform known as the kaleidoscope transform, which formalises and extends upon the concept of downsampling and concatenating an image with itself, this paper: (1) demonstrates a fundamental relationship between multiplication in modular arithmetic and downsampling; (2) provides a rigorous mathematical explanation for the fractal nature of the sampling pattern in the DFT; and (3) leverages this understanding to develop a collection of novel fractal sampling patterns for the 2D DFT with customisable properties. The ability to design tailor-made fractal sampling patterns expands the utility of the DFT in chaotic imaging and may form the basis for a bespoke chaotic sensing methodology, in which the fractal sampling matches the imaging task for improved reconstruction.

Index Terms:

Chaotic Sensing, Fractals, Kaleidoscope Transform, Sparse Image Reconstruction

I Introduction

The discrete Fourier transform (DFT) is an important and fundamental transform in many signal and imaging areas. One such area is magnetic resonance imaging (MRI), which visualises soft tissue in high resolution without the use of ionising radiation [1]; however, with long acquisition times and commensurately greater costs than other imaging modalities [2] it often goes unused.

Arguably the most impactful development towards addressing the long acquisition times of MRI is that of compressed sensing (CS) [3, 4] which, in the case of MRI, allows an image to be significantly undersampled using a random sampling pattern in the DFT of image space (often referred to as k-space) and then reconstructed with an iterative algorithm [5]; however, as in many practical applications, truly random measurements are infeasible in MRI. Fortunately, alternatives exist. Yu et al. [6] developed a hardware-friendly sampling scheme which approximates the desired properties of random sampling by using sensing matrices populated with chaotic sequences. Theirs and subsequent works have demonstrated reconstruction performance on par with or better than random sampling [7, 8, 9, 10, 11, 12].

Similarly, chaotic sensing (ChaoS) proposes the use of a deterministic fractal sampling pattern for the DFT and MRI [10]. The noise introduced by this sampling scheme is turbulent and thus image-independent, and may be easily removed through dampening and obtaining the maximum likelihood solution. Furthermore, as the fractal is constructed entirely from periodic lines in the DFT based on the discrete Radon transform [13], it may be feasibly acquired using MRI hardware [10]. ChaoS appears to be a promising alternative to CS for MRI, having been shown to outperform radial sampling and compressed sensing in structural similarity (SSIM) and visual information fidelity (VIF) when applied both to the Shepp-Logan phantom image with Gaussian noise added and to data experimentally gathered using a phantom constructed from Lego blocks in a liquid-filled tube [10]; however, although it was verified numerically that the ChaoS pattern is a fractal, it is not well understood why this fractal structure arises [10].

In this work, a novel image transform known as the kaleidoscope transform (KT) is presented, which formalises and extends upon the concept of downsampling and concatenating an image with itself. Through proving that certain scaling operations in the DFT are equivalent to KTs, the self-similar, fractal nature of the ChaoS fractal is explained. This understanding is then leveraged to develop novel, customisable fractal sampling patterns for use in ChaoS or any other methodology designed to employ the mixing of chaotic imaging information involving the DFT.

This paper is structured as follows. In the remainder of this section, the method of construction of the ChaoS fractal is provided, and the potential utility of customisable variations of this pattern is justified. Section II then introduces the KT and proves its relationship to scaling in the DFT. Finally, section III utilises the theory of the KT to explain the fractal nature of the ChaoS sampling pattern and presents novel fractal sampling patterns along with methods for their construction.

I-A Chaotic Sensing Fractal

The method for generating the fractal sampling pattern used in ChaoS on an $N\times N$ grid, where $N$ is prime, as described in [10], begins with generating the Farey sequence of order $N$ , $F_{N}$ . This is the sequence of irreducible fractions in the domain $[0,1]$ with denominator less than or equal to $N$ . This sequence may be easily generated using the mediant property of the Farey sequence. Beginning with $\frac{a_{1}}{b_{1}}=\frac{0}{1}$ and $\frac{a_{2}}{b_{2}}=\frac{1}{1}$ , the mediant fraction, $\frac{a_{3}}{b_{3}}$ , may then be generated using the rule

\displaystyle\frac{a_{3}}{b_{3}}=\left(\frac{a_{1}+a_{2}}{b_{1}+b_{2}}\right).

This may be performed recursively, stopping when $\frac{a_{3}}{b_{3}}=\frac{1}{N}$ , to generate all elements of the sequence [14].

Each of these fractions is then mapped to a grid position using the relation $\frac{a}{b}\to(b,a)$ , where $(b,a)$ represents $b$ pixels across and $a$ pixels up. It may be observed that these points represent all possible directions from the origin in the first octant of the grid. To generate the points corresponding to the rest of the plane, this octant is flipped and mirrored.

The generated points are then sorted by their Euclidean distance from the origin ( $L^{2}$ norm), $\ell_{2}=a^{2}+b^{2}$ , and those closest to the origin are selected (the exact number chosen is dependent on the Katz criterion, which will be defined shortly). Each point is then mapped to a periodic line given by all of its integer multiples computed modulo $N$ , i.e. any points on the line that would lie outside the grid are wrapped around back into it.

As mentioned, the number of periodic lines required is governed by the Katz criterion. This is given by

K=\frac{\max\left({\sum_{j=0}^{\mathcal{N}-1}|a_{j}|,\sum_{j=0}^{\mathcal{N}-1}|b_{j}|}\right)}{N},

where the minimum information required for an exact reconstruction is when $K=1$ [15]; however, $K$ may be greater than this. Using the box-counting method, it was determined that the fractal constructed with $N=257$ had a Minkowski-Bouligand dimension of 1.79; however, no attempt was made to find this dimension analytically and the reason why the self-similar, repeated copies of the central circular pattern arose from the described construction could not be explained.

I-B Fourier Signatures

Work has been conducted in the past on using the unique DFT characteristics of different image classes, known as Fourier signatures, to tailor-make MRI sampling patterns for specific tasks. Two such examples are feature-recognising MRI [16] and MRI with encoding by singular value decomposition (SVD) [17], which were found to be equivalent given the same training data [18].

Using feature-recognising MRI, $32\times 16$ pixel simulated images were able to be significantly better reconstructed than when using a direct inverse Fourier transform; however, the authors did not apply this technique to real MRI data. Furthermore, no consideration was given to the nature of the sampling determined by this method. It is foreseeable that for a given family of images it would be just as difficult to implement as the random pattern suggested by CS. Thus, although these two techniques demonstrate that there can be advantages to using sampling patterns tailored to particular image classes in MRI, neither is able to guarantee that these patterns are physically implementable on MRI hardware or that they lead to any significant acquisition acceleration.

The ChaoS fractal, on the other hand, is theoretically guaranteed to be feasibly implementable in hardware due to its construction from periodic lines in the DFT. Thus, a method to generate custom ChaoS fractals would combine the ability to tailor a sampling pattern to a specific image class, as in feature-recognising MRI, with turbulent noise of ChaoS.

II The Kaleidoscope Transform

Given a positive integer downsampling factor, $\nu$ , and integer smear factor, $\sigma$ , the $\nu,\sigma$ -KT converts a sequence, image, or array of any dimension into $\nu$ evenly spaced, downsampled copies of itself along each axis, each scaled by a factor of $\sigma$ , as presented in Fig. 1 and 2.

Refer to caption — Figure 1: $\nu,\sigma$ -KTs of a sample $512\times 512$ image. The original image (source: [19]) is presented in the $(\nu,\sigma)=(1,1)$ position.

The KT will be based on the kaleidoscope mapping and it is convenient to define two such mappings: an upper and a lower. The intuition behind these mappings is further explained in Supplementary Material A: Motivating the Kaleidoscope Mapping.

Definition 1 (Upper Kaleidoscope Mapping).

The upper $\nu,\sigma$ -kaleidoscope mapping over $\mathbb{Z}_{N}$ , where $N,\nu\in\mathbb{Z}^{+}$ and $\sigma\in\mathbb{Z}$ , written $\kappa_{U}:\mathbb{Z}_{N}\to\mathbb{Z}_{N}$ , is defined by

\displaystyle\kappa_{U}[n;\nu,\sigma]

\displaystyle=\left(\left\lceil\frac{N}{\nu}\right\rceil(n\bmod\nu)+\sigma\left\lfloor\frac{n}{\nu}\right\rfloor\right)\bmod N.

Definition 2 (Lower Kaleidoscope Mapping).

The lower $\nu,\sigma$ -kaleidoscope mapping over $\mathbb{Z}_{N}$ , where $N,\nu\in\mathbb{Z}^{+}$ and $\sigma\in\mathbb{Z}$ , written $\kappa_{L}:\mathbb{Z}_{N}\to\mathbb{Z}_{N}$ , is defined by

\displaystyle\kappa_{L}[n;\nu,\sigma]

\displaystyle=\left(\left\lfloor\frac{N}{\nu}\right\rfloor(n\bmod\nu)+\sigma\left\lfloor\frac{n}{\nu}\right\rfloor\right)\bmod N.

These may then be combined into a single, unified kaleidoscope mapping as follows, where the notation of $0^{+}$ versus $0^{-}$ is used to distinguish between cases when $\sigma=0$ .

Definition 3 (Kaleidoscope Mapping).

The $\nu,\sigma$ -kaleidoscope mapping over $\mathbb{Z}_{N}$ , where $N,\nu\in\mathbb{Z}^{+}$ and $\sigma\in\mathbb{Z}$ , written $\kappa:\mathbb{Z}_{N}\to\mathbb{Z}_{N}$ , is defined by

\displaystyle\kappa[n;\nu,\sigma]

\displaystyle=\begin{cases}\kappa_{L}[n;\nu,\sigma],&\sigma\leq 0^{-}\\ \kappa_{U}[n;\nu,\sigma],&\sigma\geq 0^{+}.\end{cases}

At first the use of this notation may appear to restrict lower and upper kaleidoscope mappings to non-positive and non-negative smear factors respectively; however, it follows simply from the definition of each mapping that smear factors that are congruent modulo $N$ give equivalent mappings. Thus, any desired smearing may be achieved with either the upper or lower mapping.

The KT moves each element of a sequence to its new position determined by the kaleidoscope mapping, summing any elements that are mapped to the same position. It is defined as follows.

Definition 4 (Kaleidoscope Transform).

The $\nu,\sigma$ -KT over $\mathbb{Z}_{N}$ , where $N,\nu\in\mathbb{Z}^{+}$ and $\sigma\in\mathbb{Z}$ , of a sequence $x[n]$ for $n\in\mathbb{Z}_{N}$ , is defined by

\displaystyle\mathcal{K}_{\nu,\sigma}\{x[n]\}[k]

\displaystyle=\sum_{n=0}^{N-1}x[n]\cdot\delta[k-\kappa[n;\nu,\sigma]]

for $k\in\mathbb{Z}_{N}$ , where $\delta$ is the unit sample function and $\kappa$ is the kaleidoscope mapping over $\mathbb{Z}_{N}$ .

This definition may easily be generalised to higher dimensions by simply performing the kaleidoscope mapping (potentially with different $\nu$ and $\sigma$ values) on the index along each dimension.

We now prove an important result that relates the kaleidoscope mapping (and hence transform) to multiplication in modular arithmetic.

Proposition 1 (Kaleidoscope-Multiplication Theorem).

If given $N,\nu\in\mathbb{Z}^{+}$ and $\sigma\in\mathbb{Z}$ such that $|\sigma|<\nu$ , we may write that $N=L\nu-\sigma$ for some $L\in\mathbb{Z^{+}}$ , then

\displaystyle\kappa[n;\nu,\sigma]

\displaystyle=Ln\bmod N,

where $\kappa[n;\nu,\sigma]$ is the $\nu,\sigma$ -kaleidoscope mapping over $\mathbb{Z}_{N}$ .

Proof.

Suppose that the premise above holds. Now, we may observe from the definition of the kaleidoscope mapping that,

\displaystyle\kappa[n;\nu,\sigma]

\displaystyle=\left(L(n\bmod\nu)+\sigma\left\lfloor\frac{n}{\nu}\right\rfloor\right)\bmod N.

Now, by the definition of the mod operator,

$\displaystyle\kappa[n;\nu,\sigma]$	$\displaystyle=\left(L\left(n-\nu\left\lfloor\frac{n}{\nu}\right\rfloor\right)+\sigma\left\lfloor\frac{n}{\nu}\right\rfloor\right)\bmod N$
	$\displaystyle=\left(Ln-\left(L\nu-\sigma\right)\left\lfloor\frac{n}{\nu}\right\rfloor\right)\bmod N$
	$\displaystyle=\left(Ln-N\left\lfloor\frac{n}{\nu}\right\rfloor\right)\bmod N$	(as $N=L\nu-\sigma$ )
	$\displaystyle=Ln\bmod N.$

∎

An important corollary of this result is that certain multiplications by an integer $L$ mod $N$ are equivalent to kaleidoscope mappings. Of particular interest are multiplications that result in unity-smear kaleidoscope mappings (those with $|\sigma|$ = 1). These occur when $L$ is a factor of $N\pm 1$ .

Multiplications by factors of $mN\pm 1$ for positive integers $m$ also approximate unity-smear kaleidoscope mappings. This is because such a multiplication mod $N$ may be expressed as the same multiplication mod $mN$ , which is a unity-smear kaleidoscope mapping of the sequence zero-padded by a factor of $m$ , followed by a reduction mod $N$ . This reduction causes each repeated pattern to fill the zero-padded gaps of its neighbours, approximating a unity-smear kaleidoscope mapping.

III Chaotic Sensing Through the Kaleidoscope

The ChaoS fractal on an $N\times N$ grid is traditionally thought of as the superposition of $R$ discrete periodic lines (where $R$ is determined by the Katz criterion), each given by $N$ multiples of an initial Farey vector; however, it may just as easily be expressed as the superposition of $\lceil N/2\rceil$ scaled images, each containing the same multiple of all $R$ Farey vectors, along with their reflections, as presented in Fig. 3. It follows simply that those scaled images containing Farey vectors multiplied by factors of $mN\pm 1$ for positive integer $m$ will approximate unity-smear KTs, containing downsampled, repeated copies of the original image (in the case when the vectors are sorted by their $L^{2}$ norm as in the classic fractal, a circle). Their superposition then gives the repeated, self-similar patterns observed in the ChaoS fractal.

This theory for the fractal nature of the ChaoS fractal predicts that the repeated pattern need not be circular, and the grid on which the fractal is constructed need not be square (only that $M\pm 1$ and $N\pm 1$ share common factors for an $M\times N$ grid). Moreover, it predicts that fractals such as these may be generated in any number of dimensions.

By using other $L^{p}$ norms to sort the Farey vectors used to generate the fractal, any superellipse may be realised as central pattern. Applying a linear transformation to each vector prior to calculating its norm allows for stretched and rotated variations of these central patterns. Furthermore, any star-shaped polygon may be realised as a central pattern by sorting the Farey vectors using what we refer to as the polygon quasi-norm. This calculates the minimum factor by which the desired polygon must be scaled for it to contain the given vector, as described in Fig. 4. Note that in practice polygon central patterns are limited to those that are 180^∘-rotationally-symmetric as the fractal is constructed from periodic lines.

Further extensions to the ChaoS fractal may be realised by diverging entirely from the paradigm of using periodic lines to construct the sampling pattern. Instead, the fractal may be explicitly constructed from the superposition of scaled copies of a desired pattern. Of particular interest are fractals constructed from feasibly implementable k-space trajectories other than lines, such as spirals. When using this method, we may also explicitly limit the scaled copies of the central pattern to those that give unity-smear KTs. Fig. 5 presents a selection of such patterns for ChaoS constructed using these methods.

IV Conclusion

In this paper the novel kaleidoscope transform was introduced, which formalises and extends upon the idea of downsampling and concatenating an image with itself. It was then demonstrated mathematically that certain scaling operations in the DFT are equivalent to KTs. This was used to explain the previous fractal nature of the ChaoS sampling pattern, and to present a variety of similar, novel patterns, along with methods for their construction for the DFT. Demonstrating that such a rich family of fractal sampling patterns exists expands upon the potential utility of the DFT for chaotic systems. In particular, it lays the foundation for a bespoke ChaoS methodology, in which the sampling pattern used is tailored to the type of object being imaged, potentially improving speed and reconstruction.

References

[1] Eric C. Ehman et al. “PET/MRI: Where might it replace PET/CT?” In Journal of magnetic resonance imaging : JMRI 46.28370695, 2017, pp. 1247–1262 URL: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5623147/
[2] Kg Hollingsworth “Reducing acquisition time in clinical MRI by data undersampling and compressed sensing reconstruction” In Physics In Medicine And Biology 60.21 IOP PUBLISHING LTD, 2015, pp. R297–R322
[3] Emmanuel Candes, Justin Romberg and Terence Tao “Stable Signal Recovery from Incomplete and Inaccurate Measurements” In Communications on Pure and Applied Mathematics 59.8 New York: John WileySons, Limited, 2006, pp. 1207–1223 URL: http://search.proquest.com/docview/214229434/
[4] D.L Donoho “Compressed sensing” In IEEE Transactions on Information Theory 52.4 IEEE, 2006, pp. 1289–1306
[5] Michael Lustig, David Donoho and John M. Pauly “Sparse MRI: The application of compressed sensing for rapid MR imaging” In Magnetic Resonance in Medicine 58.6 Hoboken: Wiley Subscription Services, Inc., A Wiley Company, 2007, pp. 1182–1195
[6] L Yu, X Zhang, W Chen and H Sun “Compressive Sensing with Chaotic Sequence” In IEEE Signal Processing Letters 17, 2010, pp. 731–734
[7] Venceslav Kafedziski and Toni Stojanovski “Compressive sampling with chaotic dynamical systems” In 2011 19th Telecommunications Forum (TELFOR) Proceedings of Papers IEEE, 2011 DOI: 10.1109/telfor.2011.6143641
[8] Li Zeng et al. “Deterministic Construction of Toeplitzed Structurally Chaotic Matrix for Compressed Sensing” In Circuits, Systems, and Signal Processing 34.3 Springer ScienceBusiness Media LLC, 2014, pp. 797–813 DOI: 10.1007/s00034-014-9873-7
[9] D. Rontani et al. “Compressive Sensing with Optical Chaos” In Scientific Reports 6.1 Springer ScienceBusiness Media LLC, 2016 DOI: 10.1038/srep35206
[10] Shekhar S Chandra et al. “Chaotic Sensing” In IEEE Transactions on Image Processing 27.12 IEEE, 2018, pp. 6079–6092
[11] Hongping Gan, Song Xiao, Tao Zhang and Feng Liu “Bipolar measurement matrix using chaotic sequence” In Communications in Nonlinear Science and Numerical Simulation 72 Elsevier BV, 2019, pp. 139–151 DOI: 10.1016/j.cnsns.2018.12.012
[12] Hongping Gan, Song Xiao and Feng Liu “Chaotic Binary Sensing Matrices” In International Journal of Bifurcation and Chaos 29.09 World Scientific Pub Co Pte Lt, 2019, pp. 1950121 DOI: 10.1142/s0218127419501219
[13] F Matúš and J Flusser “Image Representation via a Finite Radon Transform” In Pattern Analysis and Machine Intelligence, IEEE Transactions on 15.10, 1993, pp. 996–1006 DOI: 10.1109/34.254058
[14] G. Hardy and E. Wright “An Introduction to the Theory of Numbers” Oxford: Clarendon Press, 1979
[15] M.. Katz and Richard Gordon “Questions of Uniqueness and Resolution in Reconstruction from Projections” In Physics Today 32.12 American Institute of Physics, 1979, pp. 52–56
[16] Y Cao and D N Levin “Feature-recognizing MRI” In Magnetic resonance in medicine 30.3, 1993, pp. 305–317 URL: http://search.proquest.com/docview/75997850/
[17] Gary P Zientara, Lawrence P Panych and Ferenc A Jolesz “Dynamically adaptive MRI with encoding by singular value decomposition” In Magnetic resonance in medicine 32.2 Baltimore: Wiley Subscription Services, Inc., A Wiley Company, 1994, pp. 268–274
[18] Yue Cao and David N. Levin “On the Relationship Between Feature‐Recognizing MRI and MRI Encoded by Singular Value Decomposition” In Magnetic Resonance in Medicine 33.1 Baltimore: Wiley Subscription Services, Inc., A Wiley Company, 1995, pp. 140–142
[19] University of Granada Computer Vision Group “Dataset of standard 512x512 grayscale test images”, 2003 URL: http://decsai.ugr.es/cvg/CG/base.htm

Supplementary Material

IV-A Motivating the Kaleidoscope Mapping

Previously, the definitions of the kaleidoscope mapping and transform were presented with minimal justification. Here, their motivation is presented in more detail.

Suppose we have a sequence of length $N$ , $x[n]$ , which we wish to downsampled by a factor of $\nu$ (for now, we assume that $\nu|N$ ). The canonical approach is to take every $\nu^{\text{th}}$ element of $x[n]$ beginning with $x[0]$ to form the new, downsampled subsequence; however, we could just as easily have started with $x[1]$ , or $x[2]$ , or so on, all the way up to $x[\nu-1]$ . This gives us $\nu$ possible downsampled subsequences, each with different elements. We denote such a downsampled subsequence starting at the $m^{\text{th}}$ element of $x[n]$ as $y_{m}[\ell]$ . Each of these sequences has length $L=\frac{N}{\nu}$ . We may then concatenate these possible downsampled subsequences together to give the sequence $z[k]$ , also of length $N$ . An example of this process is provided in Fig. 6.

The obvious question that then arises is how we may determine which positions in $x[n]$ correspond to which positions in $z[k]$ . In other words, we seek a mapping $\kappa:n\mapsto k$ . It is relatively straightforward to see that the final index $k$ corresponding to index $\ell$ of downsampled subsequence $y_{m}[\ell]$ is given by

\displaystyle k

\displaystyle=Lm+\ell.

Noting how each downsampled sequence is constructed, it then follows that

\displaystyle k

\displaystyle=\frac{N}{\nu}(n\bmod\nu)+\left\lfloor\frac{n}{\nu}\right\rfloor.

This closely mirrors the form of the kaleidoscope mapping; however, there are some small differences. Notably, this mapping is limited to cases where $\nu|N$ . There are two approaches we may take to generalising it to any positive integer downsampling factor, $\nu$ . Firstly, we can round the length, $L$ , of each subsequence up to the nearest integer by letting $L=\left\lceil\frac{N}{\nu}\right\rceil$ . This would leave us with $\nu-1$ subsequences of equal length, and one shorter subsequence. Alternatively, we can round $L$ down to the nearest integer by letting $L=\left\lfloor\frac{N}{\nu}\right\rfloor$ . This would leave us with $\nu-1$ subsequences of equal length and one longer subsequence. As stated in the body of the paper, it is useful to consider both of these approaches, and they correspond to the upper and lower kaleidoscope mappings respectively.

We may further generalise these mappings by introducing a smear factor, $\sigma$ , which sets the spacing of each downsampled sequence when they are combined. When $\sigma=1$ , the sequences are simply concatenated, whereas other $\sigma$ values result in interleaving of the downsampled sequences, as demonstrated in Fig. 2. This generalisation is chiefly motivated by how it allows for the kaleidoscope-multiplication theorem to be stated and proved, as shown in the body of the paper.

IV-B Additional Fractal Patterns

Moving away from the context of MRI, the kaleidoscope transform may be used to generate both familiar and novel patterns, as presented in Fig. 7. Examples such as these demonstrate both the utility of the kaleidoscope transform in understanding existing self-similar patterns and its potential use in creating varied patterns for new tasks.