Cortical-Inspired Wilson–Cowan-Type Equations for Orientation-Dependent Contrast Perception Modelling

Abstract

We consider the evolution model proposed in Bertalmío (Front Comput Neurosci 8:71, 2014), Bertalmío et al. (IEEE Trans Image Process 16(4):1058–1072, 2007) to describe illusory contrast perception phenomena induced by surrounding orientations. Firstly, we highlight its analogies and differences with the widely used Wilson–Cowan equations (Wilson and Cowan in BioPhys J 12(1):1–24, 1972), mainly in terms of efficient representation properties. Then, in order to explicitly encode local directional information, we exploit the model of the primary visual cortex (V1) proposed in Citti and Sarti (J Math Imaging Vis 24(3):307–326, 2006) and largely used over the last years for several image processing problems (Duits and Franken in Q Appl Math 68(2):255–292, 2010; Prandi and Gauthier in A semidiscrete version of the Petitot model as a plausible model for anthropomorphic image reconstruction and pattern recognition. SpringerBriefs in Mathematics, Springer, Cham, 2017; Franceschiello et al. in J Math Imaging Vis 60(1):94–108, 2018). The resulting model is thus defined in the space of positions and orientation, and it is capable of describing assimilation and contrast visual bias at the same time. We report several numerical tests showing the ability of the model to reproduce, in particular, orientation-dependent phenomena such as grating induction and a modified version of the Poggendorff illusion. For this latter example, we empirically show the existence of a set of threshold parameters differentiating from inpainting to perception-type reconstructions and describing long-range connectivity between different hypercolumns in V1.

A Orientation-Dependent Model of V1

Let us denote by $R>0$ the size of the visual plane, and let $D_R\subset {{\,\mathrm{\mathbb {R}}\,}}^2$ be the disk $D_R:=\{x_1^2+x_2^2 \le R^2\}$. Fix $R>0$ such that $Q\subset D_R$. In order to exploit the properties of the roto-translation group $\mathrm{SE}(2)$ on images, we now consider them to be elements of the set:

$$\begin{aligned} \mathcal {I}= \left\{ f \in L^2({{\,\mathrm{\mathbb {R}}\,}}^2,[0,1]) \text { such that\ } {{\,\mathrm{supp}\,}}f\subset D_R\right\} . \end{aligned}$$

(24)

We remark that fixing $R>0$ is necessary, since contrast perception is strongly dependent on the scale of the features under consideration w.r.t. the visual plane.

Orientation dependence of the visual stimulus is encoded via cortical inspired techniques, following, for example, [14, 23, 27, 42, 43]. The main idea at the base of these works goes back to the 1959 paper [35] by Hubel and Wiesel (Nobel prize in 1981) who discovered the so-called hypercolumn functional architecture of the visual cortex V1.

Each neuron $\xi$ in V1 is assumed to be associated with a receptive field (RF) $\psi _\xi \in L^2({{\,\mathrm{\mathbb {R}}\,}}^2)$ such that its response under a visual stimulus $f\in \mathcal I$ is given by

$$\begin{aligned} F(\xi ) = \langle \psi _\xi , f\rangle _{L^2({{\,\mathrm{\mathbb {R}}\,}}^2)} = \int _{{{\,\mathrm{\mathbb {R}}\,}}^2} \overline{\psi _\xi ({z})} f({z})\, \mathrm{d}{z}. \end{aligned}$$

(25)

Since each neuron is sensible to a preferred position and orientation in the visual plane, we let $\xi =(x,\theta )\in \mathcal {M} = \mathbb R^2\times \mathbb P^1$. Here, $\mathbb P^1$ is the projective line that we represent as $[0,\pi ]/\sim$, with $0\sim \pi$. Moreover, in order to respect the shift-twist symmetry [17, Section  4], we will assume that the RF of different neurons are “deducible” one from the other via a linear transformation. Let us explain this in detail.

The double covering of $\mathcal {M}$ is given by the Euclidean motion group $\mathrm{SE}(2)={{\,\mathrm{\mathbb {R}}\,}}^2\rtimes \mathbb {S}^1$ that we consider endowed with its natural semi-direct product structure. That is, for $(x,\theta ),(y,\varphi )\in \mathrm{SE}(2)$, we let

$$\begin{aligned}&(x,\theta )\star (y,\varphi ) = (x+R_\theta y, \theta +\varphi ), \end{aligned}$$

(26)

$$\begin{aligned}&\text {where } R_\theta = \left( \begin{array}{ll} \cos \theta &{}\quad -\sin \theta \\ \sin \theta &{}\quad \cos \theta \end{array}\right) . \end{aligned}$$

(27)

In particular, the above operation induces an action of $\mathrm{SE}(2)$ on $\mathcal {M}$, which is thus an homogeneous space. Observe that $\mathrm{SE}(2)$ is unimodular and that its Haar measure (the left and right-invariant measure up to scalar multiples) is $\mathrm{d}x\mathrm{d}\theta$.

We now denote by $\mathcal {U}(L^2({{\,\mathrm{\mathbb {R}}\,}}^2)) \subset \mathcal {L}(L^2({{\,\mathrm{\mathbb {R}}\,}}^2))$ the space of linear unitary operators on $L^2({{\,\mathrm{\mathbb {R}}\,}}^2)$ and let ${\varPi }{:}\,\mathrm{SE}(2)\rightarrow \mathcal {U}(L^2({{\,\mathrm{\mathbb {R}}\,}}^2))$ be the quasi-regular representation of $\mathrm{SE}(2)$. That is, ${\varPi }(x,\theta )\in \mathcal {U}(L^2({{\,\mathrm{\mathbb {R}}\,}}^2))$ is the unitary operator encoding the action of the roto-translation $(x,\theta )\in \mathrm{SE}(2)$ on square-integrable functions on $\mathbb R^2$. The action of ${\varPi }(x,\theta )$ on $\psi \in L^2({{\,\mathrm{\mathbb {R}}\,}}^2)$ is

$$\begin{aligned}{}[{\varPi }(x,\theta )\psi ](y) = \psi ((x,\theta )^{-1}y) = \psi (R_{-\theta }(y-x)), \quad \forall y\in \mathbb {R}^2. \end{aligned}$$

Moreover, we let $\varLambda {:}\,\mathrm{SE}(2)\rightarrow \mathcal {U}(L^2(\mathrm{SE}(2)))$ be the left-regular representation, which acts on functions $F\in L^2(\mathrm{SE}(2))$ as

$$\begin{aligned}{}[\varLambda (x,\theta )F](y,\varphi ) = F((x,\theta )^{-1}\star (y,\varphi )) , \, \forall (y,{\varphi })\in \mathrm{SE}(2). \end{aligned}$$

(28)

Letting $L{:}\,L^2({{\,\mathrm{\mathbb {R}}\,}}^2)\rightarrow L^2(\mathcal {M})$ be the operator that transforms visual stimuli into cortical activations, one can formalise the shift-twist symmetry by requiring

$$\begin{aligned} L\circ {\varPi }(x,\theta ) = \varLambda (x,\theta )\circ L, \quad \forall (x,\theta )\in \mathrm{SE}(2). \end{aligned}$$

(29)

Under mild continuity assumption on L, it has been shown in [43] that L is then a continuous wavelet transform. That is, there exists a mother wavelet $\varPsi \in L^2({{\,\mathrm{\mathbb {R}}\,}}^2)$ satisfying ${\varPi }(x,\theta )\varPsi = {\varPi }(x,\theta +\pi )\varPsi$ for all $(x,\theta )\in \mathrm{SE}(2)$ such that

$$\begin{aligned} Lf(x,\theta ) = \langle {\varPi }(x,\theta )\varPsi , f \rangle , \quad \forall f\in L^2({{\,\mathrm{\mathbb {R}}\,}}^2), (x,\theta )\in \mathcal {M}. \end{aligned}$$

(30)

Observe that the operation ${\varPi }(x,\theta )\varPsi$ above is well defined for $(x,\theta )\in \mathcal {M}$ thanks to the assumption on $\varPsi$. By (25), the above representation of L is equivalent to the fact that the RF associated with the neuron $(x,\theta )\in \mathcal {M}$ is the roto-translation of the mother wavelet, i.e. $\psi _{(x,\theta )}={\varPi }(x,\theta )\varPsi$.

Remark 5

Letting $\varPsi ^*(x):=\overline{\varPsi (-x)}$, the above formula can be rewritten as

$$\begin{aligned} Lf(x,\theta )= & {} \int _{{{\,\mathrm{\mathbb {R}}\,}}^2} \overline{\varPsi (R_{-\theta }(y-x))}f(y)\,\mathrm{d}y \nonumber \\= & {} \big [f * (\varPsi ^*\circ R_{-\theta })\big ] (x), \quad \forall (x,\theta )\in \mathrm{SE}(2). \end{aligned}$$

(31)

where $f*g$ denotes the standard convolution on $L^2({{\,\mathrm{\mathbb {R}}\,}}^2)$.

Neurophysiological evidence shows that a good fit for the RFs is given by Gabor filters, whose Fourier transform is simply the product of a Gaussian with an oriented plane wave [25]. However, these filters are quite challenging to invert and are parametrised on a bigger space than $\mathcal M$, which takes into account also the frequency of the plane wave and not only its orientation. For this reason, in this work we chose to consider as wavelets the cake wavelets introduced in [26], see also [5]. These are obtained via a mother wavelet $\varPsi ^{\text {cake}}$ whose support in the Fourier domain is concentrated on a fixed slice, which depends on the number of orientations one aims to consider in the numerical implementation. To recover integrability properties, the Fourier transform of this mother wavelet is then smoothly cut off via a low-pass filtering, see [5, Section  2.3] for details. Observe, however, that in order to lift to $\mathcal M$ and not to $\mathrm{SE}(2)$, we consider a non-oriented version of the mother wavelet, given by $\tilde{\psi }^{\mathrm{cake}}({\omega }) + \tilde{\psi }^{\mathrm{cake}}(e^{i\pi }{\omega })$, in the notations of [5].

An important feature of cake wavelets is that, in order to recover the original image, it suffices to consider the projection operator defined by

$$\begin{aligned}&P: L^2(\mathcal M)\rightarrow L^2({{\,\mathrm{\mathbb {R}}\,}}^2), \end{aligned}$$

(32)

$$\begin{aligned}&PF(x) := \int _{\mathbb P^1} F(x,\theta )\,\mathrm{d}\theta ,\quad F\in L^2(\mathcal { M}) \end{aligned}$$

(33)

Indeed, by construction of cake wavelets, Fubini’s theorem shows that $(P\circ L)f = f$ for all $f\in \mathcal I$.