Modern Stochastics: Theory and Applications logo


  • Help
Login Register

  1. Home
  2. Issues
  3. Volume 12, Issue 1 (2025)
  4. First-order planar autoregressive model

Modern Stochastics: Theory and Applications

Submit your article Information Become a Peer-reviewer
  • Article info
  • Full article
  • Related articles
  • More
    Article info Full article Related articles

First-order planar autoregressive model
Volume 12, Issue 1 (2025), pp. 83–121
Sergiy Shklyar ORCID icon link to view author Sergiy Shklyar details  

Authors

 
Placeholder
https://doi.org/10.15559/24-VMSTA263
Pub. online: 6 August 2024      Type: Research Article      Open accessOpen Access

Received
5 February 2024
Revised
29 May 2024
Accepted
20 July 2024
Published
6 August 2024

Abstract

This paper establishes the conditions for the existence of a stationary solution to the first-order autoregressive equation on a plane as well as properties of the stationary solution. The first-order autoregressive model on a plane is defined by the equation
\[ {X_{i,j}}=a{X_{i-1,j}}+b{X_{i,j-1}}+c{X_{i-1,j-1}}+{\epsilon _{i,j}}.\]
A stationary solution X to the equation exists if and only if $(1-a-b-c)(1-a+b+c)(1+a-b+c)(1+a+b-c)\gt 0$. The stationary solution X satisfies the causality condition with respect to the white noise ϵ if and only if $1-a-b-c\gt 0$, $1-a+b+c\gt 0$, $1+a-b+c\gt 0$ and $1+a+b-c\gt 0$. A sufficient condition for X to be purely nondeterministic is provided.
An explicit expression for the autocovariance function of X on the axes is provided. With Yule–Walker equations, this facilitates the computation of the autocovariance function everywhere, at all integer points of the plane. In addition, all situations are described where different parameters determine the same autocovariance function of X.

1 Introduction

The model.  Let ${\epsilon _{i,j}}$ be a uncorrelated zero-mean equal-variance variables, $\operatorname{\mathsf{E}}{\epsilon _{i,j}}=0$ and $\operatorname{\mathsf{E}}{\epsilon _{i,j}^{2}}={\sigma _{\epsilon }^{2}}\gt 0$. Consider the equation
(1)
\[ {X_{i,j}}=a{X_{i-1,j}}+b{X_{i,j-1}}+c{X_{i-1,j-1}}+{\epsilon _{i,j}}\]
with a, b and c fixed coefficients, and X an unknown array defined on an integer lattice, for all integer indices.
Equation (1) can be used for modelling and simulating random images.
Subject of the study.  We consider the following questions about equation (1).
  • • For what coefficients a, b and c does a stationary solution of (1) exist?
  • • How to compute the autocovariance function of the stationary solution?
  • • For what coefficients a, b and c does the stationary solution satisfy the causality condition with respect to ϵ, that is, can be represented in the form
    \[ {X_{i,j}}={\sum \limits_{k=0}^{\infty }}{\sum \limits_{l=0}^{\infty }}{\psi _{k,l}}{\epsilon _{i-k,j-l}}\]
    with summable coefficients ${\psi _{k,l}}$? Actually, this question is solved in [5].
  • • How the causality condition is related to the stability of equation (1) taken nondeterministically?
  • • Is the stationary solution purely nondeterministic? This property is related to the causality condition with respect to some uncorrelated random field, called innovations: does the field of innovations exist?
  • • Which different parameters $(a,b,c,{\sigma _{\epsilon }^{2}})$ determine the same autocovariance function of the field X?
Previous study.  The review of statistical models for representing a discretely observed field on a plane, including regression models and trend-stationary models, is given in [20]. The study of planar models, their estimation and applications is presented in the monograph [7].
A theory of stationary random field on n-dimensional integer lattice was constructed by Tjøstheim [16]. For the stationary fields, he considers the property of being “purely nondeterministic,” which is related to the existence of appropriate field of innovations. He considers ARMA models, for which he provides conditions for stability of the equation and existence of a stationary solution. For AR models, he establishes Yule–Walker equations, and also consistency and asymptotic normality of the Yule–Walker estimator.
ARMA(1,1) model is studied by Basu and Reinsel in [5]. Their results on computation of the autocorrelation function can be directly applied to AR(1) model, yet in the causal case only.
Autoregressive models in a n-dimensional space had been studied before Tjøstheim. Whittle [21] studied various aspects of autoregressive model in the plane such as the spectral density of the stationary solution, nonidentifiability, and estimation of the regression coefficients. In the theoretical part of [21], the order of the model was not limited, whereas in the example a second-order model of the form
\[ {\xi _{s,t}}=\alpha {\xi _{s+1,t}}+\beta {\xi _{s-1,t}}+\gamma {\xi _{s,t+1}}+\delta {\xi _{s,t-1}}+{\epsilon _{s,t}}\]
was studied in detail.
The model equivalent to (1) is considered in [14, 19]. Two simpler models are considered in the literature. One of the models, a triangular model, is considered in the [4]. It is defined by equation ${\xi _{s,t}}=\alpha {\xi _{s+1,t}}+\beta {\xi _{s,t+1}}+{\epsilon _{s,t}}$, which, up to different notation, is also a spectral case of (1). The unit root test for this particular model is constructed in [13].
The autocovariance function of a stationary field is even-symmetric, that is ${\gamma _{X}}({h_{1}},{h_{2}})={\gamma _{X}}(-{h_{1}},-{h_{2}})$; however, ${\gamma _{X}}({h_{1}},{h_{2}})$ need not be an even function in each of its arguments. The autocovariance function of the stationary solution to (1) is even in each argument if the coefficients satisfy $c=-ab$. This special case is called a doubly geometric model. Its estimation and applications are considered in [8, 10–12]. The unit root test and estimation for parameters close to unit roots are constructed in [1, 4]. In the symmetric case, the autocovariance function is separated into two factors, one depends on ${h_{1}}$ and the other depends on ${h_{2}}$:
(2)
\[ {\gamma _{X}}({h_{1}},{h_{2}})={\gamma _{X}}({h_{1}},0){\gamma _{X}}(0,{h_{2}})/{\gamma _{X}}(0,0).\]
Random fields whose autocovariance functions satisfy separation property (2) are called linear-by-linear processes.
The property of random field of being purely nondeterministic is introduced in [16]. Further discussion of this property is found in [9] and [17].
Planar autoregressive models of order $p\times q$ (possibly of different orders along different dimensions) are considered in [6].
Methods.  We have obtained most results by manipulation with the spectral density.
Structure of the paper.  In Section 2, we provide necessary and sufficient conditions for existence of the stationary solution to equation (1). We also prove that the stationary solution must be centered and must have a spectral density. In Section 3 we study the properties of the autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$ of the stationary solution. We find ${\gamma _{X}}({h_{1}},{h_{2}})$ for ${h_{1}}=0$ and for ${h_{2}}=0$, and we prove Yule–Walker equations. This allows us to evaluate ${\gamma _{X}}({h_{1}},{h_{2}})$ everywhere. Section 4 is dedicated to the case where equation (1) is stable. In Section 5 we show nonidentifiability, there the same autocovariance function (and, in Gaussian case, the same distribution) of the random field X corresponds to two different values of the vector of parameters. In Section 6 we study Tjøstheim’s pure nondeterminism. Section 7 concludes the paper. Some auxiliary statements are presented in the appendix.

2 Existence of a stationary solution

2.1 Spectral density

Denote the shift operators
\[ {L_{x}}{X_{i,j}}={X_{i-1,j}}\hspace{1em}\text{and}\hspace{1em}{L_{y}}{X_{i,j}}={X_{i,j-1}}.\]
Equation (1) rewrites as
(3)
\[ X-a{L_{x}}X-b{L_{y}}X-c{L_{x}}{L_{y}}X=\epsilon .\]
We use the following definition of a spectral density. An integrable function ${f_{X}}({\nu _{1}},{\nu _{2}})$ is a spectral density for the wide-sense stationary field X if
(4)
\[\begin{aligned}{}{\gamma _{X}}({h_{1}},\hspace{0.1667em}{h_{2}})& =\operatorname{cov}({X_{i,j}},{X_{i-{h_{1}},j-{h_{2}}}})=\\ {} & ={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\exp (2\pi \mathrm{i}\hspace{0.1667em}({h_{1}}{\nu _{1}}+{h_{2}}{\nu _{2}})){f_{X}}({\nu _{1}},{\nu _{2}})\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}},\end{aligned}\]
where the upright i is the imaginary unit. This can be rewritten without complex numbers as
(5)
\[ {\gamma _{X}}({h_{1}},{h_{2}})={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\cos (2\pi \hspace{0.1667em}({h_{1}}{\nu _{1}}+{h_{2}}{\nu _{2}})){f_{X}}({\nu _{1}},{\nu _{2}})\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}\]
with further restraint that the spectral density should be an even function, ${f_{X}}({\nu _{1}},{\nu _{2}})={f_{X}}(-{\nu _{1}},-{\nu _{2}})$.
Due to (3), the spectral densities of the random fields ϵ and X, if they exist, are related as follows
\[ {f_{\epsilon }}({\nu _{1}},{\nu _{2}})=|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}{|^{2}}\hspace{0.1667em}{f_{X}}({\nu _{1}},{\nu _{2}}).\]
The random field ϵ is a white noise with constant spectral density ${f_{\epsilon }}({\nu _{1}},{\nu _{2}})={\sigma _{\epsilon }^{2}}$. Thus, the spectral density of the stationary field X, if it exists, is equal to
(6)
\[ {f_{X}}({\nu _{1}},{\nu _{2}})=\frac{{\sigma _{\epsilon }^{2}}}{|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}{|^{2}}}\hspace{0.1667em}.\]
Lemma 1.
Denote
\[ D=(1-a-b-c)(1-a+b+c)(1+a-b+c)(1+a+b-c).\]
Consider the denominator in (6)
\[ g({\nu _{1}},{\nu _{2}})=|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}{|^{2}},\hspace{1em}{\nu _{1}},{\nu _{2}}\in \left(-\frac{1}{2},\hspace{0.1667em}\frac{1}{2}\right].\]
It attains zero if and only if $D\le 0$. It can attain zero not more than at 2 points on the half-open box ${\left(-\frac{1}{2},\frac{1}{2}\right]^{2}}$.
Proof.
The following lines are equivalent:
(7)
\[\begin{array}{l}\displaystyle \exists {\nu _{1}}\in \left(-\frac{1}{2},\frac{1}{2}\right]\hspace{0.2778em}\exists {\nu _{2}}\in \left(-\frac{1}{2},\frac{1}{2}\right]\hspace{0.2222em}:\hspace{0.2222em}|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}{|^{2}}=0,\\ {} \displaystyle \exists {\nu _{1}}\in \left(-\frac{1}{2},\frac{1}{2}\right]\hspace{0.2778em}\exists {\nu _{2}}\in \left(-\frac{1}{2},\frac{1}{2}\right]\hspace{0.2222em}:\hspace{0.2222em}1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}={e^{2\pi \mathrm{i}{\nu _{2}}}}(b-c{e^{2\pi \mathrm{i}{\nu _{1}}}}),\end{array}\]
(8)
\[\begin{array}{l}\displaystyle \exists {\nu _{1}}\in \left(-\frac{1}{2},\frac{1}{2}\right]\hspace{0.2222em}:\hspace{0.2222em}|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}|=|(b-c{e^{2\pi \mathrm{i}{\nu _{1}}}})|,\\ {} \displaystyle \exists {\nu _{1}}\in \left(-\frac{1}{2},\frac{1}{2}\right]\hspace{0.2222em}:\hspace{0.2222em}|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}{|^{2}}=|(b-c{e^{2\pi \mathrm{i}{\nu _{1}}}}){|^{2}};\\ {} \displaystyle \exists {\nu _{1}}\in \left(-\frac{1}{2},\frac{1}{2}\right]\hspace{0.2222em}:\hspace{0.2222em}1-2a\cos (2\pi {\nu _{1}})+{a^{2}}={b^{2}}-2bc\cos (2\pi {\nu _{1}})+{c^{2}},\\ {} \displaystyle \exists {\nu _{1}}\in \left(-\frac{1}{2},\frac{1}{2}\right]\hspace{0.2222em}:\hspace{0.2222em}1+{a^{2}}-{b^{2}}-{c^{2}}=2(a-bc)\cos (2\pi {\nu _{1}}),\\ {} \displaystyle |1+{a^{2}}-{b^{2}}-{c^{2}}|\le 2\hspace{0.1667em}|a-bc|,\\ {} \displaystyle {(1+{a^{2}}-{b^{2}}-{c^{2}})^{2}}-4{(a-bc)^{2}}\le 0,\\ {} \displaystyle D\le 0.\end{array}\]
The sufficient and necessary condition for the denominator $g({\nu _{1}},{\nu _{2}})$ to attain 0 is proved.
Now treat a, b and c as fixed, and assume that $D\le 0$. Equality in (8) can be attained for not more than two ${\nu _{1}}\in \left(-\frac{1}{2},\frac{1}{2}\right]$. For each ${\nu _{1}}$, equality in (7) can be attained for no more than one ${\nu _{2}}\in \left(-\frac{1}{2},\frac{1}{2}\right]$. Thus, the function $1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}$ has not more than 2 zeros in the half-open box ${\left(-\frac{1}{2},\frac{1}{2}\right]^{2}}$. The denominator $g({\nu _{1}},{\nu _{2}})$ cannot attain zero at more than two points either.  □
It should be acknowledged that in [19] conditions for the denominator $g({\nu _{1}},{\nu _{2}})$ to attain 0 were provided, under constrains $|a|\lt 1$, $|b|\lt 1$ and $|c|\lt 1$.

2.2 Necessary and sufficient conditions for the existence of a stationary solution

Proposition 1.
Let ${\sigma ^{2}}\gt 0$, a, b and c be fixed real numbers. A collection of $\epsilon =\{{\epsilon _{i,j}},\hspace{0.2778em}i,j\in \mathbb{Z}\}$ of uncorrelated random variables with zero mean and equal variance ${\sigma ^{2}}$ and a wide-sense stationary random field $X=\{{X_{i,j}},\hspace{0.2778em}i,j\in \mathbb{Z}\}$ that satisfy equation (1) exist if and only if $D\gt 0$, where D comes from Lemma 1.
Proof.
Denote
\[ {g_{1}}({\nu _{1}},{\nu _{2}})=1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}};\]
then $g({\nu _{1}},{\nu _{2}})$ defined in Lemma 1 equals $g({\nu _{1}},{\nu _{2}})=|{g_{1}}({\nu _{1}},{\nu _{2}}){|^{2}}$. Necessity. Let random fields ϵ and X satisfy equation (1) and the other conditions of Proposition 1. Both ϵ and X are wide-sense stationary. Hence, they have spectral measures – denote them ${\lambda _{\epsilon }}$ and ${\lambda _{X}}$ – and ${\lambda _{\epsilon }}$ is a constant-weight Lebesgue measure in ${\left(-\frac{1}{2},\frac{1}{2}\right]^{2}}$. Because of (3), ${\lambda _{\epsilon }}$ has a Radon–Nikodym density $g({\nu _{1}},{\nu _{2}})$ with respect to (w.r.t.) ${\lambda _{X}}$:
\[ \frac{d{\lambda _{\epsilon }}}{d{\lambda _{X}}}({\nu _{1}},{\nu _{2}})=|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}{|^{2}}.\]
Now assume that $D\le 0$. Then, according to Lemma 1, $g({\nu _{1}},{\nu _{2}})$ attains 0 at one or two points on a half-open box ${\left(-\frac{1}{2},\frac{1}{2}\right]^{2}}$. Denote
\[\begin{array}{l}\displaystyle \begin{aligned}{}{A_{1}}& =\left\{({\nu _{1}},{\nu _{2}})\in {\left(-\frac{1}{2},\hspace{0.1667em}\frac{1}{2}\right]^{2}}:{g_{1}}({\nu _{1}},{\nu _{2}})\ne 0\right\}\\ {} & =\left\{({\nu _{1}},{\nu _{2}})\in {\left(-\frac{1}{2},\hspace{0.1667em}\frac{1}{2}\right]^{2}}:g({\nu _{1}},{\nu _{2}})\gt 0\right\},\\ {} {A_{2}}& =\left\{({\nu _{1}},{\nu _{2}})\in {\left(-\frac{1}{2},\hspace{0.1667em}\frac{1}{2}\right]^{2}}:{g_{1}}({\nu _{1}},{\nu _{2}})=0\right\}\\ {} & =\left\{({\nu _{1}},{\nu _{2}})\in {\left(-\frac{1}{2},\hspace{0.1667em}\frac{1}{2}\right]^{2}}:g({\nu _{1}},{\nu _{2}})=0\right\};\end{aligned}\\ {} \displaystyle {\lambda _{1}}(A)={\lambda _{X}}(A\cap {A_{1}}),\hspace{2em}{\lambda _{2}}(A)={\lambda _{X}}(A\cap {A_{2}})\end{array}\]
for all measurable sets $A\subset {\left(-\frac{1}{2},\hspace{0.1667em}\frac{1}{2}\right]^{2}}$.
The set ${A_{2}}$ contains one or two points. As ${\lambda _{\epsilon }}$ is absolutely continuous, ${\lambda _{\epsilon }}(A)={\lambda _{\epsilon }}(A\setminus {A_{2}})={\lambda _{\epsilon }}(A\cap {A_{1}})$. Measures ${\lambda _{1}}$ and ${\lambda _{\epsilon }}$ are absolutely continuous with respect of each other: $d{\lambda _{1}}/d{\lambda _{\epsilon }}\hspace{0.1667em}({\nu _{1}},{\nu _{2}})=g{({\mu _{1}},{\nu _{2}})^{-1}}$. Thus, ${\lambda _{1}}$ is absolutely continuous with respect to the Lebesgue measure:
\[ {\lambda _{1}}(A)={\iint _{A\cap {(-1/2,1/2]^{2}}}}\frac{{\sigma _{\epsilon }^{2}}}{g({\nu _{1}},{\nu _{2}})}\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}.\]
The measure ${\lambda _{2}}$ is concentrated at not more than 2 points. Thus, ${\lambda _{1}}$ and ${\lambda _{2}}$ are absolutely continuous and discrete components of the measure ${\lambda _{X}}$, respectively, and the nondiscrete singular component is zero.
Let $({\nu _{1}^{(0)}},{\nu _{2}^{(0)}})$ be one of the points of ${A_{2}}$. As ${g_{1}}({\nu _{1}},{\nu _{2}})$ is differentiable and $g({\nu _{1}},{\nu _{2}})\ge 0$, at the neighborhood of $({\nu _{1}^{(0)}},{\nu _{2}^{(0)}})$:
\[\begin{array}{l}\displaystyle {g_{1}}({\nu _{1}},{\nu _{2}})=o(|{\nu _{1}}-{\nu _{1}^{(0)}}|+|{\nu _{2}}-{\nu _{2}^{(0)}}|),\\ {} \displaystyle g({\nu _{1}},{\nu _{2}})=|{g_{1}}({\nu _{1}},{\nu _{2}}){|^{2}}=o({({\nu _{1}}-{\nu _{1}^{(0)}})^{2}}+{({\nu _{2}}-{\nu _{2}^{(0)}})^{2}}),\\ {} \displaystyle \frac{1}{g({\nu _{1}},{\nu _{2}})}\gt \frac{\text{const}}{{({\nu _{1}}-{\nu _{1}^{(0)}})^{2}}+{({\nu _{2}}-{\nu _{2}^{(0)}})^{2}}}\end{array}\]
for some $\text{const}\gt 0$;
\[ {\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\frac{1}{g({\nu _{1}},{\nu _{2}})}\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}=+\infty .\]
Hence
\[ {\lambda _{X}}\left({\left(-\frac{1}{2},\frac{1}{2}\right]^{2}}\right)\ge {\lambda _{1}}\left({\left(-\frac{1}{2},\frac{1}{2}\right]^{2}}\right)={\iint _{{(-1/2,1/2]^{2}}}}\frac{{\sigma ^{2}}}{g({\nu _{1}},{\nu _{2}})}\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}=\infty .\]
Thus, the spectral measure of the random field X is infinite, which is impossible. The assumption $D\le 0$ brings the contradiction. Thus, $D\gt 0$.
Sufficiency.  Let $D\gt 0$. Then ${\sigma _{\epsilon }^{2}}/g({\nu _{1}},{\nu _{2}})$ is an even integrable function on ${\left[-\frac{1}{2},\frac{1}{2}\right]^{2}}$ that attains only positive values. Then there exists a zero-mean Gaussian stationary random field X with spectral density ${\sigma _{\epsilon }^{2}}/g({\nu _{1}},{\nu _{2}})$, see Lemma 6. Define ϵ by formula (1). Then ϵ is zero-mean and stationary; ϵ has a constant spectral density ${\sigma _{\epsilon }^{2}}$. Thus, ϵ is a collection of uncorrelated random variables with zero mean and variance ${\sigma _{\epsilon }^{2}}$. Thus, the random fields ϵ and X satisfy the desired conditions.  □
Remark 1.
In the sufficiency part of Proposition 1, the probability space is specially constructed. The random fields X and ϵ constructed are jointly Gaussian and jointly stationary, which means that $\{({X_{i,j}},{\epsilon _{i,j}})\hspace{0.2778em}i,j\in \mathbb{Z}\}$ is a two-dimensional Gaussian stationary random field on a two-dimensional lattice.
Corollary 1.
Let X be a wide-sense stationary field that satisfies (1). Then X is centered and has a spectral density (6).
Proof.
Due to Proposition 1, $D\gt 0$ (where D comes from Lemma 1). Hence, $1-a-b-c\ne 0$. Taking the expectation in (1) immediately implies that the mean of the stationary random field X is zero.
Demonstration that X has a spectral density partially repeats the Necessity part of the proof of Proposition 1. In what follows, we use notations $g({\nu _{1}},{\nu _{2}})$, ${\lambda _{X}}$, ${A_{2}}$, ${\lambda _{1}}$ and ${\lambda _{2}}$ from that proof. According to Lemma 1, $D\gt 0$ implies that ${A_{2}}=\varnothing $. The discrete component ${\lambda _{2}}$ of the spectral measure ${\lambda _{X}}$ of the random filed X is concentrated on the set ${A_{2}}$, thus, ${\lambda _{2}}=0$. The nondiscrete singular component of ${\lambda _{X}}$ is also zero. The spectral measure ${\lambda _{X}}={\lambda _{1}}+{\lambda _{2}}$ is absolutely continuous, that is, the random field X has a spectral density.  □

2.3 Restatement of the existence theorem

Proposition 2.
Let ${\sigma ^{2}}\gt 0$, a, b and c be fixed real numbers, and ${\epsilon _{i,j}}$ be a collection of uncorrelated random variables with zero mean and equal variance ${\sigma ^{2}}$. A wide-sense stationary random field X that satisfies equation (1) exists if and only if $D\gt 0$, where D comes from Lemma 1.
Proof.
Necessity is logically equivalent to one in Proposition 1.
Sufficiency. Denote
(9)
\[ {\psi _{k,l}}={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\frac{\exp (-2\pi \mathrm{i}(k{\nu _{1}}+l{\nu _{2}}))}{1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}}\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}\hspace{0.1667em}.\]
The integrand is well defined due to Lemma 1. It is infinitely times differentiable periodic function. Due to Lemma 7,
(10)
\[ {\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}|{\psi _{k,l}}|\lt \infty .\]
Now prove that
(11)
\[ {X_{i,j}}={\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}{\psi _{k,l}}{\epsilon _{i-k,j-l}}\]
is a solution. The series in (11) is convergent in mean squares as well as almost surely due to Proposition 15. The resulting field X is stationary as a linear transformation of a stationary field ϵ.
By changing the indices,
\[\begin{aligned}{}{X_{i,j}}& -a{X_{i-1,j}}-b{X_{i,j-1}}-c{X_{i-1,j-1}}\\ {} & ={\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}{\psi _{k,l}}({\epsilon _{i-k,j-l}}-a{\epsilon _{i-1-k,j-l}}-b{\epsilon _{i-k,j-1-l}}-c{\epsilon _{i-1-k,j-1-l}})\\ {} & ={\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}({\psi _{k,l}}-a{\psi _{k-1,l}}-b{\psi _{k,l-1}}-c{\psi _{k-1,l-1}}){\epsilon _{i-l,j-1}}.\end{aligned}\]
This transformation of the Fourier coefficients corresponds to multiplication of the integrand by a certain function. Thus,
\[\begin{array}{l}\displaystyle \begin{aligned}{}{\psi _{k,l}}-a{\psi _{k-1,l}}& -b{\psi _{k,l-1}}-c{\psi _{k-1,l-1}}\\ {} & ={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\exp (-2\pi \mathrm{i}(k{\nu _{1}}+l{\nu _{2}}))\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}=\left\{\begin{array}{l}1\hspace{2.5pt}\text{if}\hspace{2.5pt}k=l=0,\hspace{1em}\\ {} 0\hspace{2.5pt}\text{otherwise};\hspace{1em}\end{array}\right.\end{aligned}\\ {} \displaystyle {X_{i,j}}-a{X_{i-1,j}}-b{X_{i,j-1}}-c{X_{i-1,j-1}}={\epsilon _{i,j}},\end{array}\]
and the random field X is indeed a solution to (1). Thus, X is a desired random field.  □

3 The autocovariance function

In this section the autocovariance function of the stationary process that satisfies (1) is studied. Its properties and and its values on the coordinate axes are obtained. For the causal case, which is studied in Section 4, these results follow from [5] and [2].
Assume that $\{{X_{i,j}},i,j\in \mathbb{Z}\}$ is a wide-sense stationary field that satisfies (1). According to Corollary 1, X has a spectral density, which is defined by (6). The autocovariance function can be evaluated with (4) or (5). The explicit formula is
(12)
\[\begin{aligned}{}{\gamma _{X}}({h_{1}},{h_{2}})& ={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\frac{{\sigma _{\epsilon }^{2}}\exp (2\pi \mathrm{i}({h_{1}}{\nu _{1}}+{h_{2}}{\nu _{2}}))}{|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}{|^{2}}}\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}\\ {} & ={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\frac{{\sigma _{\epsilon }^{2}}\exp (2\pi \mathrm{i}({h_{1}}{\nu _{1}}+{h_{2}}{\nu _{2}}))}{g({\nu _{1}},{\nu _{2}})}\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}},\end{aligned}\]
with $g({\nu _{1}},{\nu _{2}})$ defined in Lemma 1.
For fixed ${h_{1}}$ and ${h_{2}}$, it is possible to compute the integral in (12) as a function of a, b and c.

3.1 The autocovariance function for ${h_{1}}=0$

The denominator $g({\nu _{1}},{\nu _{2}})$ in (6) is equal to
\[\begin{aligned}{}g({\nu _{1}},{\nu _{2}})& =1+{a^{2}}+{b^{2}}+{c^{2}}+2(ac-b)\cos (2\pi {\nu _{2}})\\ {} & \hspace{1em}+2((bc-a)+(ab-c)\cos (2\pi {\nu _{2}}))\cos (2\pi {\nu _{1}})\\ {} & \hspace{1em}+2(ab+c)\sin (2\pi {\nu _{1}})\sin (2\pi {\nu _{2}})\\ {} & =A+B\cos (2\pi {\nu _{1}})+C\sin (2\pi {\nu _{1}}),\end{aligned}\]
with
(13)
\[ \begin{array}{c}\displaystyle A=1+{a^{2}}+{b^{2}}+{c^{2}}+2(ac-b)\cos (2\pi {\nu _{2}}),\\ {} \displaystyle B=2((bc-a)+(ab-c)\cos (2\pi {\nu _{2}})),\hspace{2em}C=2(ab+c)\sin (2\pi {\nu _{2}}).\end{array}\]
We are going to use Lemma 8 to calculate the integral ${\textstyle\int _{-1/2}^{1/2}}g{({\nu _{1}},{\nu _{2}})^{-1}}\hspace{0.1667em}d{\nu _{1}}$. We have to verify the conditions $A\gt 0$ and ${A^{2}}-{B^{2}}-{C^{2}}\gt 0$:
(14)
\[ {(1+{a^{2}}+{b^{2}}+{c^{2}})^{2}}-4{(ac-b)^{2}}=({(a-c)^{2}}+{(b+1)^{2}})({(a+c)^{2}}+{(b-1)^{2}})\ge 0.\]
The zero is attained if either $a=c$ and $b=-1$, or $a+c=0$ and $b=1$. In both cases $D=0$, while $D\gt 0$ according to Proposition 1 (D is defined in Lemma 1). Thus the inequality in (14) is strict,
\[\begin{array}{l}\displaystyle {(1+{a^{2}}+{b^{2}}+{c^{2}})^{2}}\gt 4{(ac-b)^{2}},\\ {} \displaystyle 1+{a^{2}}+{b^{2}}+{c^{2}}\gt |2(ac-b)|\ge -2(ac-b)\cos (2\pi {\nu _{2}}),\\ {} \displaystyle A=1+{a^{2}}+{b^{2}}+{c^{2}}+2(ac-b)\cos (2\pi {\nu _{2}})\gt 0.\end{array}\]
With some trigonometric transformations,
\[ {A^{2}}-{B^{2}}-{C^{2}}={(1-{a^{2}}+{b^{2}}-{c^{2}}-2(ac+b)\cos (2\pi {\nu _{2}}))^{2}}.\]
Again, with some transformations and having in mind that $D\gt 0$,
\[\begin{array}{l}\displaystyle {(1-{a^{2}}+{b^{2}}-{c^{2}})^{2}}-4{(ac+b)^{2}}=D\gt 0,\\ {} \displaystyle |1-{a^{2}}+{b^{2}}-{c^{2}}|\gt |2(ac+b)|\ge |2(ac+b)\cos (2\pi {\nu _{2}})|,\\ {} \displaystyle 1-{a^{2}}+{b^{2}}-{c^{2}}\ne 2(ac+b)\cos (2\pi {\nu _{2}}),\\ {} \displaystyle {A^{2}}-{B^{2}}-{C^{2}}={(1-{a^{2}}+{b^{2}}-{c^{2}}-2(ac+b)\cos (2\pi {\nu _{2}}))^{2}}\gt 0.\end{array}\]
According to Lemma 8, equation (42),
\[ {\int _{-1/2}^{1/2}}\frac{d{\nu _{1}}}{g({\nu _{1}},{\nu _{2}})}=\frac{1}{\sqrt{{A^{2}}-{B^{2}}-{C^{2}}}}=\frac{1}{|1-{a^{2}}+{b^{2}}-{c^{2}}-2(ac+b)\cos (2\pi {\nu _{2}})|}.\]
Denote
\[ {g_{2}}({\nu _{2}})=1-{a^{2}}+{b^{2}}-{c^{2}}-2(ac+b)\cos (2\pi {\nu _{2}})\]
and compute ${\textstyle\int _{-1/2}^{1/2}}{e^{2\pi \mathrm{i}{\nu _{2}}{h_{2}}}}\hspace{0.1667em}|{g_{2}}({\nu _{2}}){|^{-1}}\hspace{0.1667em}d{\nu _{2}}$ using formula (42). The expression ${g_{2}}({\nu _{2}})$ does not attain 0; hence, either ${g_{2}}({\nu _{2}})\gt 0$ for all ${\nu _{2}}$, or ${g_{2}}({\nu _{2}})\lt 0$ for all ${\nu _{2}}$. Then
\[ |{g_{2}}({\nu _{2}})|={A_{2}}+{B_{2}}\cos (2\pi {\nu _{2}}),\]
where ${A_{2}}=1-{a^{2}}+{b^{2}}-{c^{2}}$ and ${B_{2}}=-2(ac+b)$ if ${g_{2}}({\nu _{2}})\gt 0$ for all ${\nu _{2}}$; otherwise, ${A_{2}}=-(1-{a^{2}}+{b^{2}}-{c^{2}})$ and ${B_{2}}=2(ac+b)$ if ${g_{2}}({\nu _{2}})\lt 0$ for all ${\nu _{2}}$. In both cases, ${A_{2}}=|{g_{2}}(\pi /2)|\gt 0$ and
\[ {A_{2}^{2}}-{B_{2}^{2}}={(1-{a^{2}}+{b^{2}}-{c^{2}})^{2}}-4{(ac+b)^{2}}=D\gt 0.\]
According to Lemma 8, equation (43),
\[ {\int _{-1/2}^{1/2}}\frac{\exp (2\pi \mathrm{i}{h_{2}}{\nu _{2}})}{|{g_{2}}({\nu _{2}})|}\hspace{0.1667em}d{\nu _{2}}=\frac{{\beta ^{|{h_{2}}|}}}{\sqrt{{A_{2}^{2}}-{B_{2}^{2}}}}=\frac{1}{\sqrt{D}}\]
where
\[ \beta =\frac{{B_{2}}}{{A_{2}}+\sqrt{{A_{2}^{2}}-{B_{2}^{2}}}}=\frac{{B_{2}}}{{A_{2}}+\sqrt{D}}.\]
The value of β is
\[\begin{aligned}{}\beta & =0\hspace{1em}\text{if}\hspace{2.5pt}ac+b=0,\\ {} \beta & =\frac{1-{a^{2}}+{b^{2}}-{c^{2}}-\sqrt{D}}{2(ac+b)}\hspace{1em}\text{if}\hspace{2.5pt}ac+b\ne 0\hspace{2.5pt}\text{and}\hspace{2.5pt}{g_{2}}({\nu _{2}})\gt 0\hspace{2.5pt}\text{for all}\hspace{2.5pt}{\nu _{2}},\\ {} \beta & =\frac{1-{a^{2}}+{b^{2}}-{c^{2}}+\sqrt{D}}{2(ac+b)}\hspace{1em}\text{if}\hspace{2.5pt}ac+b\ne 0\hspace{2.5pt}\text{and}\hspace{2.5pt}{g_{2}}({\nu _{2}})\lt 0\hspace{2.5pt}\text{for all}\hspace{2.5pt}{\nu _{2}};\end{aligned}\]
the formula that if valid in all cases is
(15)
\[ \beta =\frac{2(ac+b)}{1-{a^{2}}+{b^{2}}-{c^{2}}+\operatorname{sign}(1-{a^{2}}+{b^{2}}-{c^{2}})\sqrt{D}}.\]
Finally,
(16)
\[\begin{aligned}{}{\gamma _{X}}(0,{h_{2}})& ={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\frac{{\sigma _{\epsilon }^{2}}\exp (2\pi \mathrm{i}{h_{2}}{\nu _{2}})}{g({\nu _{1}},{\nu _{2}})}\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}\\ {} & ={\int _{-1/2}^{1/2}}\frac{{\sigma _{\epsilon }^{2}}\exp (2\pi \mathrm{i}{h_{2}}{\nu _{2}})}{|{g_{2}}({\nu _{2}})|}\hspace{0.1667em}d{\nu _{2}}=\frac{{\beta ^{|{h_{2}}|}}{\sigma _{\epsilon }^{2}}}{\sqrt{D}}.\end{aligned}\]
In particular,
(17)
\[ \operatorname{var}{X_{i,j}}={\gamma _{X}}(0,0)=\frac{{\sigma _{\epsilon }^{2}}}{\sqrt{D}}.\]
Due to symmetry,
(18)
\[ {\gamma _{X}}({h_{1}},0)=\frac{{\alpha ^{|{h_{1}}|}}{\sigma _{\epsilon }^{2}}}{\sqrt{D}},\]
where
(19)
\[ \alpha =\frac{2(a+bc)}{1+{a^{2}}-{b^{2}}-{c^{2}}+\operatorname{sign}(1+{a^{2}}-{b^{2}}-{c^{2}})\sqrt{D}}.\]
Notice that $|\alpha |\lt 1$ and $|\beta |\lt 1$, whence
(20)
\[ |{\gamma _{X}}(1,0)|\lt {\gamma _{X}}(0,0)\hspace{1em}\text{and}\hspace{1em}|{\gamma _{X}}(0,1)|\lt {\gamma _{X}}(0,0).\]

3.2 Yule–Walker equations

We formally write down Yule–Walker equations for the autocovariance function of X. These equations are particular cases of ones given in [16, Section 6]:
(21)
\[\begin{aligned}{}{\gamma _{X}}(0,\hspace{0.1667em}0)& =a{\gamma _{X}}(1,\hspace{0.1667em}0)+b{\gamma _{X}}(0,\hspace{0.1667em}1)+c{\gamma _{X}}(1,\hspace{0.1667em}1)+{\sigma _{\epsilon }^{2}},\end{aligned}\]
(22)
\[\begin{aligned}{}{\gamma _{X}}({h_{1}},\hspace{0.1667em}{h_{2}})& =a{\gamma _{X}}({h_{1}}-1,\hspace{0.1667em}{h_{2}})+b{\gamma _{X}}({h_{1}},\hspace{0.1667em}{h_{2}}-1)+c{\gamma _{X}}({h_{1}}-1,\hspace{0.1667em}{h_{2}}-1).\end{aligned}\]
In Lemma 2 we obtain conditions for (22) to hold true. As (21) requires specific conditions on coefficients a, b and c, we postpone the consideration of (21) to Section 4.
Lemma 2.
Let X be a stationary field satisfying equation (1), and let ${\gamma _{X}}({h_{1}},{h_{2}})$ be the covariance function of the process X. Then equality (22) holds true if any of the following conditions hold true:
  • (1) ${h_{1}}\ge 1$, $(1-a-b-c)(1+a-b+c)\gt 0$ and $(1-a+b+c)(1+a+b-c)\gt 0$,
  • (2) ${h_{1}}\le 0$, $(1-a-b-c)(1+a-b+c)\lt 0$ and $(1-a+b+c)(1+a+b-c)\lt 0$,
  • (3) ${h_{2}}\ge 1$, $(1-a-b-c)(1-a+b+c)\gt 0$ and $(1+a-b+c)(1+a+b-c)\gt 0$, or
  • (4) ${h_{2}}\le 0$, $(1-a-b-c)(1-a+b+c)\lt 0$ and $(1+a-b+c)(1+a+b-c)\lt 0$.
Proof.
According to Corollary 1, the stationary field X has a spectral density, that is (4) with ${f_{X}}({h_{1}},{h_{2}})$ defined in (6). Hence, according to rules how the Fourier transform changes when the function is linearly transformed,
(23)
\[\begin{array}{cc}& \displaystyle {\gamma _{X}}({h_{1}},\hspace{0.1667em}{h_{2}})-a{\gamma _{X}}({h_{1}}-1,\hspace{0.1667em}{h_{2}})-b\gamma ({h_{1}},\hspace{0.1667em}{h_{2}}-1)-c\gamma ({h_{1}}-1,\hspace{0.1667em}{h_{2}}-1)\\ {} & \displaystyle \begin{aligned}{}& ={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}{e^{2\pi \mathrm{i}({h_{1}}{\nu _{1}}+{h_{2}}{\nu _{2}})}}(1-a{e^{-2\pi \mathrm{i}{\nu _{1}}}}-b{e^{-2\pi \mathrm{i}{\nu _{2}}}}-c{e^{-2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}})\\ {} & \hspace{2em}\hspace{2em}\hspace{2em}\times {f_{X}}({h_{1}},{h_{2}})\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}\\ {} & ={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\frac{\exp (2\pi \mathrm{i}({h_{1}}{\nu _{1}}+{h_{2}}{\nu _{2}})){\sigma _{\epsilon }^{2}}}{1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}}\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}.\end{aligned}\end{array}\]
Case (1): ${h_{1}}\ge 1$, $(1-a-b-c)(1+a-b+c)\gt 0$ and $(1-a+b+c)(1+a+b-c)\gt 0$. In this case,
\[\begin{array}{l}\displaystyle {(1-b)^{2}}\gt {(a+c)^{2}}\hspace{1em}\text{and}\hspace{1em}{(1+b)^{2}}\gt {(a-c)^{2}},\\ {} \displaystyle 1-{a^{2}}+{b^{2}}-{c^{2}}\gt 2b+2ac\hspace{1em}\text{and}\hspace{1em}1-{a^{2}}+{b^{2}}-{c^{2}}\gt -2b-2ac,\\ {} \displaystyle 1-{a^{2}}+{b^{2}}-{c^{2}}\gt 2\hspace{0.1667em}|b+ac|.\end{array}\]
Then, for every ${\nu _{2}}\in [-1/2,\hspace{0.2222em}1/2]$, since $|\cos (2\pi {\nu _{2}})|\le 1$,
(24)
\[\begin{array}{l}\displaystyle 1-{a^{2}}+{b^{2}}-{c^{2}}\gt 2(b+ac)\cos (2\pi {\nu _{2}}),\\ {} \displaystyle 1-2b\cos (2\pi {\nu _{2}})+{b^{2}}\gt {a^{2}}+2ac\cos (2\pi {\nu _{2}})+{c^{2}},\\ {} \displaystyle |1-b{e^{2\pi \mathrm{i}{\nu _{2}}}}{|^{2}}\gt |a+c{e^{2\pi \mathrm{i}{\nu _{2}}}}{|^{2}},\\ {} \displaystyle |1-b{e^{2\pi \mathrm{i}{\nu _{2}}}}|\gt |a+c{e^{2\pi \mathrm{i}{\nu _{2}}}}|.\end{array}\]
According to Lemma 9,
(25)
\[ {\int _{-1/2}^{1/2}}\frac{\exp (2\pi \mathrm{i}{h_{1}}{\nu _{1}})}{1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}}\hspace{0.1667em}d{\nu _{1}}=0,\]
which, with (23), implies (22).
Case (2): ${h_{1}}\le 0$, $(1-a-b-c)(1+a-b+c)\lt 0$ and $(1-a+b+c)(1+a+b-c)\lt 0$. In this case,
\[\begin{array}{l}\displaystyle {(1-b)^{2}}\lt {(a+c)^{2}}\hspace{1em}\text{and}\hspace{1em}{(1+b)^{2}}\lt {(a-c)^{2}},\\ {} \displaystyle -2b-2ac\lt {a^{2}}-{b^{2}}+{c^{2}}-1\hspace{1em}\text{and}\hspace{1em}2b+2ac\lt {a^{2}}-{b^{2}}+{c^{2}}-1,\\ {} \displaystyle 2\hspace{0.1667em}|b+ac|\lt {a^{2}}-{b^{2}}+{c^{2}}-1.\end{array}\]
Then, for every ${\nu _{2}}\in [-1/2,\hspace{0.2222em}1/2]$, since $|\cos (2\pi {\nu _{2}})|\le 1$,
\[\begin{array}{l}\displaystyle -2(b+ac)\cos (2\pi {\nu _{2}})\le 2\hspace{0.1667em}|b+ac|\lt {a^{2}}-{b^{2}}+{c^{2}}-1,\\ {} \displaystyle 1-2b\cos (2\pi {\nu _{2}})+{b^{2}}\lt {a^{2}}+2ac\cos (2\pi {\nu _{2}})+{c^{2}},\\ {} \displaystyle |1-b{e^{2\pi \mathrm{i}{\nu _{2}}}}{|^{2}}\lt |a+c{e^{2\pi \mathrm{i}{\nu _{2}}}}{|^{2}},\\ {} \displaystyle |1-b{e^{2\pi \mathrm{i}{\nu _{2}}}}|\lt |a+c{e^{2\pi \mathrm{i}{\nu _{2}}}}|.\end{array}\]
Again, according to another case of Lemma 9, (25) still holds true. With (23), it implies (22).
Cases (3) and (4) are symmetric to cases (1) and (2), respectively.  □

3.3 Uniqueness of the solution

We prove that under Sufficient conditions of Proposition 2, the stationary solution to (1) is unique.
Proposition 3.
Let ${\sigma ^{2}}\gt 0$, a, b and c be fixed real numbers, and ${\epsilon _{i,j}}$ be a collection of uncorrelated random variables with zero mean and equal variance ${\sigma ^{2}}$. If $D\gt 0$, then equation (1) has only one stationary solution X, namely the one defined by (11) with coefficients ${\psi _{k,l}}$ defined in (9). The coefficients ${\psi _{k,l}}$ satisfy (10); hence, the double series in (11) converges in least squares and almost surely.
Proof.
It remains to prove the uniqueness of the solution, as other assertions follow from the proof of Proposition 2.
Let X be a solution to (1), which can be rewritten as (3). Let $\mathcal{H}$ be the Hilbert space spanned by ${L_{x}^{k}}{L_{y}^{l}}X$, $k,l=\dots ,-1,0,1,\dots $, with scalar product
\[ \langle {Y^{(1)}},{Y^{(2)}}\rangle =\operatorname{\mathsf{E}}{Y_{0,0}^{(1)}}{Y_{0,0}^{(2)}}.\]
Formally, $\mathcal{H}$ can be defined as the closure of the set of random fields of form $\{{\textstyle\sum _{k=-n}^{n}}{\textstyle\sum _{l=-n}^{n}}{\phi _{k,l}}{X_{i-k,j-l}},\hspace{0.2778em}i,j=\dots ,-1,0,1,\dots \}$ in a Banach space of bounded-second-moment random fields with norm
(26)
\[ \| Y\| =\underset{i,j}{\sup }{(\operatorname{\mathsf{E}}{Y_{i,j}^{2}})^{1/2}}.\]
All random fields within $\mathcal{H}$ are jointly stationary, which imply that the scalar product indeed corresponds to the norm of the Banach space:
\[ \langle Y,Y\rangle =\| Y{\| ^{2}}\hspace{1em}\text{for all}\hspace{1em}Y\in \mathcal{H}.\]
Since (3), $\epsilon \in \mathcal{H}$. Let us construct a similar space for ϵ. Let
\[\begin{aligned}{}& {\mathcal{H}_{1}}=\left\{{\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}{\theta _{k,l}}{L_{x}^{k}}{L_{y}^{k}}\epsilon :{\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}{\theta _{k,l}^{2}}\lt \infty \right\}\\ {} & =\left\{\left\{{\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}{\theta _{k,l}}{\epsilon _{i-k,j-l}},\hspace{0.1667em}i,j=\dots ,-1,0,1,\dots \hspace{-0.1667em}\right\}:{\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}{\theta _{k,l}^{2}}\lt \infty \hspace{-0.1667em}\right\}.\end{aligned}\]
The series is convergent, see Proposition 16.
The random fields ${\sigma _{\epsilon }^{-1}}{L_{x}^{k}}{L_{y}^{k}}\epsilon $ make an orthonormal basis of the subspace ${\mathcal{H}_{1}}$. The orthogonal projector onto ${\mathcal{H}_{1}}$ can be defined as follows:
\[\begin{array}{l}\displaystyle {P_{1}}Y=\frac{1}{{\sigma _{\epsilon }^{2}}}{\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}(\operatorname{\mathsf{E}}{Y_{k,l}}{\epsilon _{0,0}}){L_{x}^{k}}{L_{y}^{l}}\epsilon ,\\ {} \displaystyle {({P_{1}}Y)_{i,j}}=\frac{1}{{\sigma _{\epsilon }^{2}}}{\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}\operatorname{\mathsf{E}}({Y_{k,l}}{\epsilon _{0,0}}){\epsilon _{i-k,j-l}}\end{array}\]
for all $Y\in \mathcal{H}$.
It is possible to verify that
\[ {P_{1}}{L_{x}}={L_{x}}{P_{1}},\hspace{2em}{P_{1}}{L_{y}}={L_{y}}{P_{1}},\hspace{2em}{P_{1}}\epsilon =\epsilon .\]
Applying the operator ${P_{1}}$ to both sides of (3), we get that the random field ${P_{1}}X$ is also a stationary solution to (1). Both random fields must be centered, and variances of both the fields satisfy (17):
\[ \| X{\| ^{2}}=\| {P_{1}}X\| =\frac{{\sigma _{\epsilon }^{2}}}{\sqrt{D}}.\]
This implies that $\| X-{P_{1}}X\| =0$, which means $X={P_{1}}X$ almost surely.
Combining (10) and (23), we get
(27)
\[\begin{array}{l}\displaystyle {\gamma _{X}}(-k,-l)-a{\gamma _{X}}(-k-1,-l)-b{\gamma _{X}}(-k,-l-1)-c{\gamma _{X}}(-k-1,-l-1)={\sigma _{\epsilon }^{2}}{\psi _{k,l}},\\ {} \displaystyle \operatorname{\mathsf{E}}{X_{0,0}}{X_{k,l}}-a\operatorname{\mathsf{E}}{X_{-1,0}}{X_{k,l}}-b\operatorname{\mathsf{E}}{X_{0,-1}}{X_{k,l}}-c\operatorname{\mathsf{E}}{X_{-1,-1}}{X_{k,l}}={\sigma _{\epsilon }^{2}}{\psi _{k,l}},\\ {} \displaystyle \operatorname{\mathsf{E}}({X_{0,0}}-a{X_{-1,0}}-b{X_{0,-1}}-c{X_{-1,-1}}){X_{k,l}}={\sigma _{\epsilon }^{2}}{\psi _{k,l}},\\ {} \displaystyle \operatorname{\mathsf{E}}{\epsilon _{0,0}}{X_{k,l}}={\sigma _{\epsilon }^{2}}{\psi _{k,l}},\\ {} \displaystyle {X_{i,j}}=P{X_{i,j}}=\frac{1}{{\sigma _{\epsilon }^{2}}}{\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}(\operatorname{\mathsf{E}}{X_{k,l}}{\epsilon _{0,0}}){\epsilon _{i-k,j-l}}={\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}{\psi _{k,l}}{\epsilon _{i-k,j-l}}.\end{array}\]
Thus, the solution to (1) must satisfy (11) almost surely, whence the uniqueness follows.  □

4 Stability and causality

4.1 Causality

We generalize the causality condition for random fields on a two-dimensional lattice.
Definition 1.
Let ${\epsilon _{i,j}}$ (i and j integers) be random variables of zero mean and constant, finite, nonzero variance. A random field $\{{X_{i,j}},\hspace{0.2778em}i,j=\dots ,-1,0,1,\dots \}$ is said to be causal with respect to (w.r.t.) the white noise ϵ if there exists an array of coefficients $\{{\psi _{k,l}},\hspace{0.2778em}k,l=0,1,2,\dots \}$ such that
(28)
\[ {\sum \limits_{k=0}^{\infty }}{\sum \limits_{l=0}^{\infty }}|{\psi _{k,l}}|\lt \infty \]
and the process X allows the representation
(29)
\[ {X_{i,j}}={\sum \limits_{k=0}^{\infty }}{\sum \limits_{l=0}^{\infty }}{\psi _{k,l}}{\epsilon _{i-k,j-l}}.\]
A causal random field must be wide-sense stationary.
Proposition 4.
Let ${\epsilon _{i,j}}$ be a uncorrelated zero-mean equal-variance variables, $\operatorname{\mathsf{E}}{\epsilon _{i,j}}=0$ and $\operatorname{\mathsf{E}}{\epsilon _{i,j}^{2}}={\sigma _{\epsilon }^{2}}$, $0\lt {\sigma _{\epsilon }^{2}}\lt \infty $. Let X be a stationary field that satisfies (1). The random field X is causal w.r.t. the white noise ϵ if and only if these four inequalities hold true:
\[\begin{aligned}{}1-a-b-c& \gt 0,\hspace{2em}1-a+b+c\gt 0,\\ {} 1+a-b+c& \gt 0,\hspace{2em}1+a+b-c\gt 0.\end{aligned}\]
If the field X is causal, then the coefficients in representation (29) equal
(30)
\[ {\psi _{k,l}}={\sum \limits_{m=0}^{\min (k,l)}}\left(\genfrac{}{}{0.0pt}{}{k}{m}\right)\left(\genfrac{}{}{0.0pt}{}{l}{m}\right){a^{k-m}}{b^{l-m}}{(ab+c)^{m}}.\]
Proof.
Sufficiency. Method 1. Denote
(31)
\[ {f_{1}}=1-a-b-c,\hspace{1em}{f_{2}}=1-a+b+c,\hspace{1em}{f_{3}}=1+a-b+c,\hspace{1em}{f_{4}}=1+a+b-c.\]
Obviously, ${f_{1}}+{f_{2}}+{f_{3}}+{f_{4}}=4\gt 0$.
The conditions are ${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$ and ${f_{4}}\gt 0$; assume they hold true. Verify the conditions in [5, Proposition 1]:
\[ 1-|\alpha |=\min \hspace{-0.1667em}\left(\frac{{f_{1}}+{f_{2}}}{2},\hspace{0.2222em}\frac{{f_{3}}+{f_{4}}}{2}\right)\gt 0,\]
whence $|a|\lt 1$. Similarly, $|b|\lt 1$ and $|c|\lt 1$. Also,
\[\begin{array}{l}\displaystyle {(1+{a^{2}}-{b^{2}}-{c^{2}})^{2}}-4{(a+bc)^{2}}={f_{1}}{f_{2}}{f_{3}}{f_{4}}\gt 0,\\ {} \displaystyle 1-{b^{2}}-|a+bc|=\frac{{f_{1}}{f_{4}}+{f_{2}}{f_{3}}+2\min ({f_{1}}{f_{2}},\hspace{0.1667em}{f_{3}}{f_{4}})}{4}\gt 0.\end{array}\]
Thus, conditions of [5, Proposition 1] are satisfied. Hence, the polynomial function $\Phi ({z_{1}},{z_{2}})=1-a{z_{1}}-b{z_{2}}-c{z_{1}}{z_{2}}$ does not attain zero on the set $\{({z_{1}},{z_{2}})\in {\mathbb{C}^{2}}:|{z_{1}}|\le 1\text{and}|{z_{2}}|\le 1\}$. According to [16, Theorem 5.1], the function $\Psi {({z_{1}},{z_{2}})^{-1}}$ allows the Taylor expansion
\[ \Phi ({z_{1}},{z_{2}})={\sum \limits_{k=0}^{\infty }}{\sum \limits_{l=0}^{\infty }}{\psi _{k,l}}{z_{1}^{k}}{z_{2}^{l}}\]
with summable coefficients. The process X defined in (29) satisfies (1). With uniqueness stated in Proposition 3, this completes the proof.
Method 2. The stationary solution X to equation (1) admits a representation (11), with the coefficients ${\psi _{k,l}}$ defined in (9); the coefficients satisfy (10). With application of Lemma 9 as it has been done in the proof of Lemma 2, it is possible to demonstrate that ${\psi _{k,l}}=0$ if either $k\lt 0$ or $l\lt 0$. Thus, (28) and (29) hold true.
Necessity. Method 1. Let the process X admit representation (29) with condition (28) being satisfied. Now use formalism, which originates in [5, Section 5]. Denote a bivariate power series $\Psi ({z_{1}},{z_{2}})={\textstyle\sum _{k=0}^{\infty }}{\textstyle\sum _{l=0}^{\infty }}{\psi _{k,l}}{z_{1}^{k}}{z_{2}^{l}}$. The function ${\Psi _{(}}{z_{1}},{z_{2}})$ is bounded on the set $\{({z_{1}},{z_{2}})\in {\mathbb{C}^{2}}:|{z_{1}}|\le 1\text{and}|{z_{2}}|\le 1\}$. Then $X=\Psi ({L_{x}},{L_{y}})\epsilon $, and equation $\Phi ({L_{x}},{L_{y}})\Psi ({L_{x}},{L_{y}})\epsilon =\epsilon $ holds true, with polynomial $\Phi ({z_{1}},{z_{2}})=1-a{z_{1}}-b{z_{2}}-c{z_{1}}{z_{2}}$ defined in the Sufficiency part. This is possible only if $\Phi ({L_{x}},{L_{y}})\Psi ({L_{x}},{L_{y}})$ is the identity operator, and $\Phi ({z_{1}},{z_{2}})\Psi ({z_{1}},{z_{2}})=1$. (The operators ${L_{x}}$ and ${L_{y}}$ commute and, as they are considered linear operators acting on the Banach space of bounded-second-moment random fields with norm (26), their norm equals 1. Hence, the substitution of the operators ${L_{x}}$ and ${L_{y}}$ into a formal power series with summable coefficients makes sense.) Hence, $\Phi ({z_{1}},{z_{2}})$ is nonzero for any complex numbers ${z_{1}}$ and ${z_{2}}$ such that $|{z_{1}}|\le 1$ and $|{z_{1}}|\le 1$. Necessary and sufficient conditions stated in [5, Proposition 1] must be satisfied, in particular, $|a|\lt 1$, $|b|\lt 1$, $|c|\lt 1$ and ${f_{1}}{f_{2}}{f_{3}}{f_{4}}={(1+{a^{2}}-{b^{2}}-{c^{2}})^{2}}-4{(a+bc)^{2}}\gt 0$. It could be verified that if any two of the factors ${f_{1}}$, ${f_{2}}$, ${f_{3}}$ or ${f_{4}}$ are negative, then $\max (|a|,|b|,|c|)\gt 1$, which is false. Thus, conditions in [5] imply that ${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$ and ${f_{4}}\gt 0$.
Method 2. According to Lemma 1, $D={f_{1}}{f_{2}}{f_{3}}{f_{4}}\gt 0$. Thus, the sum and the product of the factors ${f_{1}}$, ${f_{2}}$, ${f_{3}}$ and ${f_{4}}$ are positive. This is possible in two cases: (1) all ${f_{1}}$, ${f_{2}}$, ${f_{3}}$ and ${f_{4}}$ are positive; or (2) two of the factors ${f_{1}}$, ${f_{2}}$, ${f_{3}}$ and ${f_{4}}$ are positive, and the other two are negative. Case (1) is what is to be proved. In case (2) coefficients defined in (9) are zero, ${\psi _{k,l}}=0$ for all $k\ge 0$ and $l=0$, as follows from the proof of Lemma 2. Thus, all coefficients ${\psi _{k,l}}$ in the representation (11) are zero, which is absurd.
The expression for ${\psi _{k,l}}$. Again, if ${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$ and ${f_{4}}\gt 0$, in the proof of Lemma 2, case (1), inequality (24) was demonstrated. Let k and l be nonnegative integers. By (9) and Lemmas 9 and 10,
\[\begin{aligned}{}{\psi _{k,l}}& ={\int _{-1/2}^{1/2}}{e^{-2\pi \mathrm{i}l{\nu _{2}}}}\bigg({\int _{-1/2}^{1/2}}\frac{\exp (-2\pi \mathrm{i}k{\nu _{1}})}{1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}}\hspace{0.1667em}d{\nu _{1}}\bigg)\hspace{0.1667em}d{\nu _{2}}\\ {} & ={\int _{-1/2}^{1/2}}\frac{{e^{-2\pi \mathrm{i}l{\nu _{2}}}}{(a+c{e^{2\pi \mathrm{i}{\nu _{2}}}})^{k}}}{{(1-b{e^{2\pi \mathrm{i}{\nu _{2}}}})^{k+1}}}\hspace{0.1667em}d{\nu _{2}}\\ {} & ={\sum \limits_{n=0}^{\min (k,l)}}\left(\genfrac{}{}{0.0pt}{}{k}{n}\right)\left(\genfrac{}{}{0.0pt}{}{l}{n}\right){a^{k-n}}{b^{l-n}}{(ab+c)^{n}}.\end{aligned}\]
 □
Remark 2.
The equivalent simplification of causality or stability conditions stated in [5] was done in [19, 2, 3]. The set of points whose coordinates $(a,b,c)$ satisfy the causality condition is found in all these papers.
Notice that occasionally the set of necessary and sufficient conditions contains a redundant condition, which follows from other conditions. For example, in [5, Proposition] the inequality $1-{b^{2}}\gt |a+bc|$ follows from $|a|\lt 1$, $|b|\lt 1$, $|c|\lt 1$ and ${(1+{a^{2}}-{b^{2}}-{c^{2}})^{2}}-4{(a+bc)^{2}}\gt 0$. In [2] inequalities $|a|\lt 1$, $|b|\lt 1$ and $|c|\lt 1$ follow from $a-b-c\lt 1$, $-a+b-c\lt 1$, $-a-b+c\lt 1$ and $a+b+c\lt 1$.
The conditions for causality can be restated as follows. Let X be a stationary random field that satisfies (1). Then X is causal with respect to the white noise ϵ if and only if $|a|\lt 1$, $|b|\lt 1$ and $|c|\lt 1$. The trick is that $D\gt 0$ (the necessary and sufficient condition for existence of such a field X, see Proposition 1) is assumed.
Remark 3.
Refs. [5, eq. (2.1)] and [2, eq. (1.4)–(1.5)] provide simpler and explicit expressions for the coefficients ${\psi _{k,l}}$:
\[\begin{aligned}{}{\psi _{k,l}}& ={\sum \limits_{m=0}^{\min (k,l)}}\frac{(k+l-m)!}{(k-m)!\hspace{0.1667em}(k-m)!\hspace{0.1667em}m!}{a^{k-m}}{b^{l-m}}{c^{m}},\\ {} {\psi _{k,l}}& =\left(\genfrac{}{}{0.0pt}{}{k+l}{l}\right){a^{k}}{b^{l}}{_{2}}\hspace{-0.1667em}{F_{1}}\hspace{-0.1667em}\left(-k,-l;-k-l;-\frac{c}{ab}\right)\hspace{1em}\text{if}\hspace{2.5pt}ab\ne 0,\end{aligned}\]
where ${_{2}}\hspace{-0.1667em}{F_{1}}$ is the Gauss hypergeometric function.

4.2 Autocovariance function and Yule–Walker equations under the causality condition

The following Proposition can be proved in the same way as a similar proposition for time series in one dimension, so the proof is skipped.
Proposition 5.
Let coefficients a, b and c be such that ${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$ and ${f_{4}}\gt 0$ for ${f_{1}},\dots ,{f_{4}}$ defined in (31). Let X be a stationary field and ϵ be a collection of zero-mean equal-variance random variables, $\operatorname{var}{\epsilon _{i,j}}={\sigma _{\epsilon }^{2}}$, that satisfy (1). Then $\operatorname{cov}({X_{i,j}},{\epsilon _{k,l}})=0$ if $k\gt i$ or $l\gt j$, and the following Yule–Walker equations hold true:
\[\begin{aligned}{}{\gamma _{X}}(0,0)& =a{\gamma _{X}}(1,0)+b{\gamma _{X}}(0,1)+c{\gamma _{X}}(1,1)+{\sigma _{\epsilon }^{2}},\\ {} {\gamma _{X}}({h_{1}},{h_{2}})& =a{\gamma _{X}}({h_{1}}-1,{h_{2}})+b{\gamma _{X}}({h_{1}},{h_{2}}-1)+c{\gamma _{X}}({h_{1}}-1,{h_{2}}-1)\end{aligned}\]
if $\max ({h_{1}},{h_{2}})\gt 0$.
Here $\operatorname{cov}({X_{i,j}},{\epsilon _{k,l}})=0$ (for $k\gt i$ or $l\gt j$) immediately follows from causality, and (22) (for $\max ({h_{1}},{h_{2}})\ge 0$) immediately follows from Lemma 2. The latter equation was also proved in [5].
Proposition 6.
Let coefficients a, b and c be such that ${f_{1}}\gt 0,\hspace{0.2222em}\dots ,\hspace{0.2222em}{f_{4}}\gt 0$, with ${f_{1}},\dots ,{f_{4}}$ defined in (31). Let X be a stationary random field, and ϵ be a collection of zero-mean random variables with equal variance $\operatorname{\mathsf{E}}{\epsilon _{i,j}^{2}}$ that satisfy (1). Then random fields X and ϵ are jointly stationary with cross-autocovariance function
\[\begin{array}{l}\displaystyle \operatorname{cov}({X_{i+{h_{1}},j+{h_{2}}}},{\epsilon _{i,j}})=0\hspace{1em}\textit{if}\hspace{2.5pt}{h_{1}}\lt 0\hspace{2.5pt}\textit{or}\hspace{2.5pt}{h_{2}}\lt 0,\\ {} \displaystyle \operatorname{cov}({X_{i+{h_{1}},j+{h_{2}}}},{\epsilon _{i,j}})={\sum \limits_{k=0}^{\min ({h_{1}},{h_{2}})}}\left(\genfrac{}{}{0.0pt}{}{{h_{1}}}{k}\right)\left(\genfrac{}{}{0.0pt}{}{{h_{2}}}{k}\right){a^{{h_{1}}-k}}{b^{{h_{2}}-k}}{(ab+c)^{k}}{\sigma _{\epsilon }^{2}}\end{array}\]
if ${h_{1}}\ge 0$ and ${h_{2}}\ge 0$.
Proof.
The joint stationarity follows from (3) and the stationarity of X (this was already used in the proof of Proposition 3). The formula for the autocovariance function follows from (27) and Proposition 4.  □
The following proposition provides a simple explicit expression for the autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$ for arguments ${h_{1}}$ and ${h_{2}}$ of opposite sign for the causal case. These formulas are proved in [5, Proposition 1] for a more general model (namely, for ARMA(1,1) model).
Proposition 7.
Under conditions of Proposition 6, the autocovariance function of the field X satisfies
(32)
\[ {\gamma _{X}}({h_{1}},{h_{2}})=\frac{{\alpha ^{|{h_{1}}|}}{\beta ^{|{h_{2}}|}}{\sigma _{\epsilon }^{2}}}{\sqrt{D}}\hspace{1em}\textit{if}\hspace{2.5pt}{h_{1}}{h_{2}}\le 0\textit{,}\]
where D is defined in Lemma 1, and
\[ \alpha =\frac{2(a+bc)}{1+{a^{2}}-{b^{2}}-{c^{2}}+\sqrt{D}},\hspace{2em}\beta =\frac{2(ac+b)}{1-{a^{2}}+{b^{2}}-{c^{2}}+\sqrt{D}}.\]
The proof is by induction, by using results of Section 3.1 for the base case and Yule–Walker equations for the induction step. As that adds nothing to [5], we skip the detailed proof.
Thus, if ${h_{1}}{h_{2}}\le 0$, then
\[\begin{array}{l}\displaystyle {\gamma _{X}}({h_{1}},{h_{2}}){\gamma _{X}}(0,0)={\gamma _{X}}({h_{1}},0){\gamma _{X}}(0,{h_{2}}),\\ {} \displaystyle \operatorname{corr}({X_{i,j}},{X_{i-{h_{1}},j-{h_{2}}}})=\operatorname{corr}({X_{i,j}},{X_{i-{h_{1}},j}})\operatorname{corr}({X_{i,j}},{X_{i,j-{h_{2}}}}).\end{array}\]
This result can be generalized to autoregressive models of higher order, see [6].
Under conditions of Proposition 7,
(33)
\[ {\gamma _{X}}(1,-1){\gamma _{X}}(0,0)={\gamma _{X}}(1,0){\gamma _{X}}(0,1).\]
The results of Section 4.2 allow to compute the autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$ recursively. The Taylor expansion of the autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$, with respect to the parameter c, is provided in [19].

4.3 Stability

Lemma 3.
Let a, b and c be such that ${f_{1}},\dots ,{f_{4}}$ defined in (31) are positive real numbers. Let ${\psi _{k,l}}$ be defined in (30) if $k\ge 0$ and $l\ge 0$, and ${\psi _{k,l}}=0$ if $k\lt 0$ or $l\lt 0$. Then ${\psi _{k,l}}$ satisfy this equation:
(34)
\[ {\psi _{k,l}}-a{\psi _{k-1,l}}-b{\psi _{k,l-1}}-c{\psi _{k-1,l-1}}=\left\{\begin{array}{l}1\hspace{1em}\textit{if}\hspace{2.5pt}k=l=0,\hspace{1em}\\ {} 0\hspace{1em}\textit{otherwise}.\hspace{1em}\end{array}\right.\]
Proof.
Under conditions ${f_{1}}\gt 0,\hspace{0.1667em}\dots ,\hspace{0.1667em}{f_{4}}\gt 0$, the coefficients ${\psi _{k,l}}$ satisfy (9) for all integer k and l. Respective transforms of the Fourier coefficients and the integrand in (9) yield
\[\begin{array}{c}\displaystyle {\psi _{k,l}}-a{\psi _{k-1,l}}-b{\psi _{k,l-1}}-c{\psi _{k-1,l-1}}\\ {} \displaystyle ={\int _{-1/2}^{1/2}}{\int _{-1/2}^{1/2}}\exp (-2\pi \mathrm{i}(k{\nu _{1}}+l{\nu _{2}}))\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}=\left\{\begin{array}{l}1\hspace{1em}\text{if}\hspace{2.5pt}k=l=0,\hspace{1em}\\ {} 0\hspace{1em}\text{otherwise}.\hspace{1em}\end{array}\right.\end{array}\]
 □
Remark 4.
In Lemma 3 the condition $D\gt 0$ is unnecessary. Indeed, condition $D\gt 0$ holds true for small a, b and c. The coefficients ${\psi _{k,l}}$ are polynomials in a, b and c; thus, the left-hand side of (34) is also a polynomial in a, b and c, and the fact that (34) holds true for small a, b and c implies that (34) holds true for all real a, b and c.
Proposition 8.
Let a, b and c be such that ${f_{1}},\dots ,{f_{4}}$ defined in (31) are positive real numbers. Consider the recurrence equation
(35)
\[ {x_{i,j}}=a{x_{i-1,j}}+b{x_{i,j-1}}+c{x_{i-1,j-1}}+{v_{i,j}}\]
for all positive integers i and j, where
  • • ${v_{i,j}}$ are known variables, $i\gt 0$ and $j\gt 0$;
  • • ${x_{0,0}}$, ${x_{i,0}}$ and ${x_{0,j}}$ are preset, so setting their values makes initial/boundary conditions;
  • • ${x_{i,j}}$ ($i\gt 0$ and $j\gt 0$) are unknown variables.
Then the explicit formula for the solution is the following:
(36)
\[\begin{aligned}{}{x_{i,j}}& ={\psi _{i,j}}{x_{0,0}}+{\sum \limits_{k=0}^{i-1}}{\psi _{k,j}}\hspace{0.1667em}({x_{i-k,0}}-a{x_{i-k-1,0}})\\ {} & \hspace{1em}+{\sum \limits_{l=0}^{j-1}}{\psi _{i,l}}\hspace{0.1667em}({x_{0,j-l}}-b{x_{0,j-l-1}})+{\sum \limits_{k=0}^{i-1}}{\sum \limits_{l=0}^{j-1}}{\psi _{k,l}}{v_{i-k,j-l}},\end{aligned}\]
with ${\psi _{k,l}}$ defined in (30).
If, in addition, the sequences of ${x_{i,0}}$ and ${x_{0,j}}$ and the array of ${v_{i,j}}$ are bounded, then the solution x is also bounded.
Proof.
First, show that under convention that “an empty sum is zero,” namely ${\textstyle\sum _{i=0}^{-1}}\dots =0$ and ${\textstyle\sum _{j=0}^{-1}}\dots =0$, the right side of (36) satisfies the boundary conditions. For one or both indices being zero, the coefficients ${\psi _{k,l}}$ equal
\[ {\psi _{0,0}}=1,\hspace{2em}{\psi _{k,0}}={a^{k}},\hspace{2em}{\psi _{0,l}}={b^{l}}\]
for all integer $k\ge 0$ and $l\ge 0$. Hence,
\[\begin{aligned}{}{\psi _{0,0}}{x_{0,0}}& ={x_{i,j}},\\ {} {\psi _{0,j}}{x_{0,0}}& +{\sum \limits_{l=0}^{j-1}}{\psi _{0,l}}({x_{0,j-l}}-b{x_{0,j-l-1}})={b^{l}}{x_{0,0}}+{\sum \limits_{l=0}^{j-1}}{b^{l}}({x_{0,j-k}}-b{x_{0,j-l-1}})\\ {} & ={b^{l}}{x_{0,0}}+{\sum \limits_{l=0}^{j-1}}{b^{l}}{x_{0,j-l}}-{\sum \limits_{l=0}^{j-1}}{b^{l+1}}{x_{0,j-l-1}}\\ {} & ={b^{l}}{x_{0,0}}+{\sum \limits_{l=0}^{j-1}}{b^{l}}{x_{0,j-l}}-{\sum \limits_{l=1}^{j}}{b^{l}}{x_{0,j-l}}={b^{l}}{x_{0,0}}+{x_{0,j}}-{b^{l}}{x_{0,0}}={x_{0,j}}\end{aligned}\]
for integer $j\gt 0$, and due to symmetry for all integers $i\gt 0$
\[ {\psi _{i,0}}{x_{0,0}}+{\sum \limits_{k=0}^{i-1}}{\psi _{k,0}}\hspace{0.1667em}({x_{i-k,0}}-a{x_{i-k-1,0}})={x_{i,0}}.\]
Known variables ${v_{i,j}}$ are defined for all integer $i\gt 0$ and $j\gt 0$. Extend their domain for all integer $i\ge 0$ and $j\ge 0$ by assigning
\[ {v_{0,0}}=0;\hspace{2em}{v_{i,0}}={x_{i,0}}-a{x_{i-1,0}},\hspace{1em}i\gt 0;\hspace{2em}{v_{0,j}}={x_{0,j}}-b{x_{0,j-1}},\hspace{1em}j\gt 0.\]
Then the right-hand side of (36) is
\[ {x_{i,j}^{\mathrm{RHS}}}={\sum \limits_{k=0}^{i}}{\sum \limits_{l=0}^{j}}{\psi _{k,l}}{v_{i-k,j-l}},\hspace{1em}i\ge 0,\hspace{2.5pt}j\ge 0.\]
In particular, ${x_{i,j}^{\mathrm{RHS}}}={x_{i,j}}$ if either $i=0$ or $j=0$.
Now prove that ${x^{\mathrm{RHS}}}$ satisfies (35). To that end, let ${\psi _{k,l}}=0$ if $k\lt 0$ or $l\lt 0$. Then, for i and j positive integers,
\[\begin{array}{c}\displaystyle {x_{i,j}^{\mathrm{RHS}}}-a{x_{i-1,j}^{\mathrm{RHS}}}-b{x_{i,j-1}^{\mathrm{RHS}}}-c{x_{i-1,j-1}^{\mathrm{RHS}}}\\ {} \displaystyle \begin{aligned}{}& ={\sum \limits_{k=0}^{i}}{\sum \limits_{l=0}^{j}}{\psi _{k,l}}{v_{i-k,j-l}}-a{\sum \limits_{k=0}^{i-1}}{\sum \limits_{l=0}^{j}}{\psi _{k,l}}{v_{i-1-k,j-l}}\\ {} & \hspace{1em}-b{\sum \limits_{k=0}^{i}}{\sum \limits_{l=0}^{j-1}}{\psi _{k,l}}{v_{i-k,j-1-l}}-c{\sum \limits_{k=0}^{i-1}}{\sum \limits_{l=0}^{j-1}}{\psi _{k,l}}{v_{i-1-k,j-1-l}}\\ {} & ={\sum \limits_{k=0}^{i}}{\sum \limits_{l=0}^{j}}{\psi _{k,l}}{v_{i-k,j-l}}-a{\sum \limits_{k=1}^{i}}{\sum \limits_{l=0}^{j}}{\psi _{k-1,l}}{v_{i-k,j-l}}\\ {} & \hspace{1em}-b{\sum \limits_{k=0}^{i}}{\sum \limits_{l=1}^{j}}{\psi _{k,l-1}}{v_{i-k,j-l}}-c{\sum \limits_{k=1}^{i}}{\sum \limits_{l=1}^{j}}{\psi _{k-1,l-1}}{v_{i-k,j-l}}\\ {} & ={\sum \limits_{k=0}^{i}}{\sum \limits_{l=0}^{j}}({\psi _{k,l}}-a{\psi _{k-1,l}}-b{\psi _{k,l-1}}-c{\psi _{k-1,l-1}}){v_{i-k,j-l}}={v_{i,j}}\hspace{0.1667em};\end{aligned}\end{array}\]
here we used Lemma 3.
Finally, equality ${x_{i,j}}={x_{i,j}^{\mathrm{RHS}}}$ for all integer $i\ge 0$ and $j\ge 0$ can be proved by induction.
If initial/boundary values and ${v_{i,j}}$ are bounded,
\[ \underset{i\gt 0}{\sup }|{x_{i,0}}|\lt \infty ,\hspace{1em}\underset{j\gt 0}{\sup }|{x_{0,j}}|\lt \infty ,\hspace{1em}\text{and}\hspace{1em}\underset{i,j\gt 00}{\sup }|{v_{i,j}}|\lt \infty ,\]
then
\[\begin{aligned}{}|{x_{i,j}}|& =\left|{\sum \limits_{k=0}^{i}}{\sum \limits_{k=0}^{j}}{\psi _{k,l}}{v_{i-k,j-l}}\right|\le {\sum \limits_{k=0}^{i}}{\sum \limits_{j=0}^{l}}|{\psi _{i,j}}|\underset{0\le r\le 1}{\max }\underset{0\le s\le 1}{\max }|{v_{r,s}}|\\ {} & \le {\sum \limits_{k=0}^{i}}{\sum \limits_{j=0}^{l}}|{\psi _{i,j}}|\times \begin{aligned}{}\max \Big(& |{x_{0,0}}|,\hspace{0.2778em}(1+|a|)\underset{1\le r\le i}{\max }|{x_{r,0}}|,\\ {} & (1+|b|)\underset{1\le s\le i}{\max }|{x_{0,s}}|,\hspace{0.2778em}\underset{1\le r\le i}{\max }\underset{1\le s\le i}{\max }|{v_{i,j}}|\Big)\end{aligned}\\ {} & \le {\sum \limits_{k=0}^{\infty }}{\sum \limits_{j=0}^{\infty }}|{\psi _{i,j}}|\times \begin{aligned}{}\max \Big(& |{x_{0,0}}|,\hspace{0.2778em}(1+|a|)\underset{r\ge 1}{\sup }|{x_{r,0}}|,\\ {} & (1+|b|)\underset{s\ge 1}{\sup }|{x_{0,s}}|,\hspace{0.2778em}\underset{r,s\ge 1}{\sup }|{v_{i,j}}|\Big)\lt \infty .\end{aligned}\end{aligned}\]
Thus, the solution x is bounded. Here we used (28).  □

5 Symmetry, nonidentifiability, and special cases

In this section we consider the question how parameters change when the random field X is flipped, and whether two or more sets of parameters correspond to the same distribution (or the same autocovariance function) of X. Here are four parameters: three coefficients a, b, and c and the variance of the error term ${\sigma _{\epsilon }^{2}}$. Usually, two sets of parameters $(a,b,c,{\sigma _{\epsilon }^{2}})$ correspond to the same autocovariance function. In special cases, up to four sets of parameters can correspond to the same autocovariance function. We start with special cases.

5.1 Special cases

The purpose of this section is to illustrate the use of results of Sections 2–4. The symmetric case was studied in [10, 5]. In that case, the random field X is also called a linear-by-linear process or a doubly-geometric process. For example, the formula for the autocorrelation function in the symmetric case similar to (37) in Proposition 10 was given in [10, Section 3].

5.1.1 Symmetric case: $ab+c=0$

Proposition 9.
Let a stationary field X and a collection of uncorrelated zero-mean equal-variance random variables ϵ be a solution to (1). The autocovariance function ${\gamma _{X}}$ is even with respect to each variable, ${\gamma _{X}}({h_{1}},{h_{2}})={\gamma _{X}}(-{h_{1}},{h_{2}})={\gamma _{X}}({h_{1}},-{h_{2}})$ if and only if $ab+c=0$.
Proof.
Equation (4) can be rewritten as
\[ {\gamma _{X}}({h_{1}},\hspace{0.1667em}{h_{2}})={\int _{-1/2}^{1/2}}{e^{2\pi \mathrm{i}{h_{2}}{\nu _{2}}}}{\int _{-1/2}^{1/2}}{e^{2\pi \mathrm{i}{h_{1}}{\nu _{1}}}}{f_{X}}({\nu _{1}},{\nu _{2}})\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}.\]
The even symmetry with respect to ${h_{2}}$, ${\gamma _{X}}({h_{1}},{h_{2}})={\gamma _{X}}({h_{1}},-{h_{2}})$, is equivalent to the inner integral being a real function,
\[\begin{aligned}{}& {\int _{-1/2}^{1/2}}{e^{2\pi \mathrm{i}{h_{1}}{\nu _{1}}}}{f_{X}}({\nu _{1}},{\nu _{2}})\hspace{0.1667em}d{\nu _{1}}\\ {} & ={\int _{-1/2}^{1/2}}\frac{\exp (2\pi \hspace{0.1667em}\mathrm{i}\hspace{0.1667em}{h_{1}}{\nu _{1}}){\sigma _{\epsilon }^{2}}}{|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}-b{e^{2\pi \mathrm{i}{\nu _{2}}}}-c{e^{2\pi \mathrm{i}({\nu _{1}}+{\nu _{2}})}}{|^{2}}}\hspace{0.1667em}d{\nu _{1}}\\ {} & ={\int _{-1/2}^{1/2}}\frac{\exp (2\pi \hspace{0.1667em}\mathrm{i}\hspace{0.1667em}{h_{1}}{\nu _{1}}){\sigma _{\epsilon }^{2}}}{A+B\cos (2\pi {\nu _{1}})+C\sin (2\pi {\nu _{2}})}\in \mathbb{R}\hspace{2.5pt}\text{for almost all}\hspace{2.5pt}{\nu _{2}}\in \left[-\frac{1}{2},\frac{1}{2}\right],\end{aligned}\]
with A, B and C defined in (13) in Section 3.1. In turn, this is equivalent to the function Fourier-transformed being an even function,
\[ \frac{{\sigma _{\epsilon }^{2}}}{A+B\cos (2\pi {\nu _{1}})+C\sin (2\pi {\nu _{1}})}=\frac{{\sigma _{\epsilon }^{2}}}{A+B\cos (-2\pi {\nu _{1}})+C\sin (-2\pi {\nu _{1}})}\]
for almost all ${\nu _{1}}$ and ${\nu _{2}}\in \left[-\frac{1}{2},\frac{1}{2}\right]$, which holds true if and only if $C=2(ab+c)\sin (2\pi {\nu _{2}})=0$ for all ${\nu _{2}}$, which is equivalent to $ab+c=0$.  □
Thus, we consider the case $c=-ab$. The following proposition establishes conditions on a and b.
Proposition 10.
Let $c=-ab$. The stationary field X and a collection of uncorrelated zero-mean variables with equal variance ${\sigma _{\epsilon }^{2}}$ which satisfy (1) exist if and only if $|a|\ne 1$ and $|b|\ne 1$. In that case, the autocovariance function of X is equal
(37)
\[ {\gamma _{X}}({h_{1}},{h_{2}})=\frac{{a^{\pm {h_{1}}}}{b^{\pm {h_{2}}}}{\sigma _{\epsilon }^{2}}}{|1-{a^{2}}|\times |1-{b^{2}}|},\]
where the signs “±” are chosen so that $|{a^{\pm {h_{2}}}}|\lt 1$ and $|{b^{\pm {h_{2}}}}|\lt 1$.
The field X is causal with respect to ϵ if and only if $|a|\lt 1$ and $|b|\lt 1$.
Proof.
Apply Proposition 1. Since $c=-ab$,
\[\begin{aligned}{}D& =(1-a-b+ab)(1-a+b-ab)(1+a-b-ab)(1+a+b+ab)=\\ {} & ={(1-a)^{2}}{(1-b)^{2}}{(1+b)^{2}}{(1+a)^{2}}.\end{aligned}\]
Thus, $D\ge 0$, and the necessary and sufficient condition $D\gt 0$ for the desired process to exist is satisfied if and only if $(1-a)(1-b)(1+b)(1+a)\ne 0$, which is equivalent to $|a|\ne 1$ and $|b|\ne 1$.
The spectral density (6) splits into two factors:
\[ {f_{X}}({\nu _{1}},{\nu _{2}})=\frac{1}{|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}{|^{2}}}\times \frac{{\sigma _{\epsilon }^{2}}}{|1-b{e^{2\pi \mathrm{i}{\nu _{2}}}}{|^{2}}}.\]
Thus, the autocovariance function is
\[ {\gamma _{X}}({h_{1}},\hspace{0.1667em}{h_{2}})={\int _{-1/2}^{1/2}}\frac{\exp (2\pi {h_{1}}{\nu _{1}})}{|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}{|^{2}}}\hspace{0.1667em}d{\nu _{1}}\times {\int _{-1/2}^{1/2}}\frac{\exp (2\pi {h_{2}}{\nu _{2}})}{|1-b{e^{2\pi \mathrm{i}{\nu _{2}}}}{|^{2}}}\hspace{0.1667em}d{\nu _{2}}\hspace{0.1667em}{\sigma _{\epsilon }^{2}}.\]
With Lemma 8,
\[\begin{aligned}{}{\int _{-1/2}^{1/2}}\frac{\exp (2\pi {h_{1}}{\nu _{1}})}{|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}{|^{2}}}\hspace{0.1667em}d{\nu _{1}}& ={\int _{-1/2}^{1/2}}\frac{\exp (2\pi {h_{1}}{\nu _{1}})}{1+{a^{2}}-2a\cos (2\pi {\nu _{1}})}\hspace{0.1667em}d{\nu _{1}}\\ {} & =\frac{{\alpha ^{|{h_{1}}|}}}{\sqrt{{(1+{a^{2}})^{2}}-4{a^{2}}}}=\frac{{\alpha ^{|{h_{1}}|}}}{|1-{a^{2}}|}\end{aligned}\]
with
\[ \alpha =\frac{2a}{1+{a^{2}}+|1-{a^{2}}|}=\left\{\begin{array}{l@{\hskip10.0pt}l}a\hspace{1em}& \text{if}\hspace{2.5pt}|a|\lt 1,\\ {} {a^{-1}}\hspace{1em}& \text{if}\hspace{2.5pt}|a|\gt 1.\end{array}\right.\]
Thus,
\[ {\int _{-1/2}^{1/2}}\frac{\exp (2\pi {h_{1}}{\nu _{1}})}{|1-a{e^{2\pi \mathrm{i}{\nu _{1}}}}{|^{2}}}\hspace{0.1667em}d{\nu _{1}}=\frac{{a^{\pm {h_{1}}}}}{|1-{a^{2}}|}.\]
Similarly,
\[ {\int _{-1/2}^{1/2}}\frac{\exp (2\pi {h_{2}}{\nu _{2}})}{|1-b{e^{2\pi \mathrm{i}{\nu _{2}}}}{|^{2}}}\hspace{0.1667em}d{\nu _{2}}=\frac{{b^{\pm {h_{2}}}}}{|1-{b^{2}}|},\]
and (37) holds true.
The causality conditions from Proposition 4 rewrite as follows:
\[\begin{aligned}{}(1-a)(1-b)& \gt 0,\hspace{2em}(1-a)(1+b)\gt 0,\\ {} (1+a)(1-b)& \gt 0,\hspace{2em}(1+a)(1+b)\gt 0.\end{aligned}\]
They mean that $1-a$, $1-b$, $1+a$ and $1+b$ should be (nonzero and) of the same sign. As sum of these factors $(1-a)+(1-b)+(1+a)+(1+b)=4$ is positive, the causality condition is equivalent to
\[ 1-a\gt 0,\hspace{1em}1-b\gt 0,\hspace{1em}1+a\gt 0,\hspace{1em}\text{and}\hspace{1em}1+b\gt 0,\]
which in turn is equivalent to $|a|\lt 1$ and $|b|\lt 1$.  □

5.1.2 Case $a=-bc$: uncorrelated observations along transect line

Another special case is where the coefficients in (1) satisfy relation $a=-bc$.
Proposition 11.
Let $a=-bc$. The stationary field X and a collection of uncorrelated zero-mean variables with equal variance ${\sigma _{\epsilon }^{2}}$ which satisfy (1) exist if and only if $|b|\ne 1$ and $|c|\ne 1$. The autocovariance function of X satisfies relations
\[ {\gamma _{X}}({h_{1}},0)=0,\hspace{1em}{h_{1}}\ne 0\hspace{0.1667em},\hspace{2em}{\gamma _{X}}(0,{h_{2}})=\frac{{b^{\pm {h_{2}}}}{\sigma _{\epsilon }^{2}}}{|1-{b^{2}}|\times |1-{c^{2}}|}\hspace{0.1667em},\]
where the sign “±” is chosen so that $|{b^{\pm {h_{2}}}}|\lt 1$.
The field X is causal with respect to ϵ if and only if $|b|\lt 1$ and $|c|\lt 1$.
Proof.
Since $a=-ac$,
\[\begin{aligned}{}D& =(1+bc-b-c)(1+bc+b+c)(1-bc-b+c)(1-bc+b-c)=\\ {} & ={(1-b)^{2}}{(1-c)^{2}}{(1+b)^{2}}{(1+c)^{2}}.\end{aligned}\]
Thus, $D\ge 0$, and the necessary and sufficient condition of existence of a stationary solution, namely $D\gt 0$, is satisfied if and only if $(1-b)(1-c)(1+b)(1+c)\ne 0$, which is equivalent to $|b|\ne 1$ and $|c|\ne 1$.
For the autocovariance function, we use results of Section 3.1. With notations (15) and (19),
\[\begin{array}{l}\displaystyle \sqrt{D}=|1-{b^{2}}|\times |1-{c^{2}}|,\hspace{2em}\alpha =0,\\ {} \displaystyle \begin{aligned}{}\beta & =\frac{2(-{b^{2}}c+b)}{1-{b^{2}}{c^{2}}+{b^{2}}-{c^{2}}+\operatorname{sign}(1-{b^{2}}{c^{2}}+{b^{2}}-{c^{2}})\sqrt{D}}\\ {} & =\frac{2(1-{c^{2}})b}{(1+{b^{2}})(1-{c^{2}})+\operatorname{sign}(1-{c^{2}})|1-{b^{2}}|\times |1-{c^{2}}|}\\ {} & =\frac{2(1-{c^{2}})b}{(1+{b^{2}})(1-{c^{2}})+|1-{b^{2}}|\hspace{0.1667em}(1-{c^{2}})}\\ {} & =\frac{2b}{(1+{b^{2}})+|1-{b^{2}}|}=\left\{\begin{array}{l@{\hskip10.0pt}l}b\hspace{1em}& \text{if}\hspace{2.5pt}|b|\lt 1\text{,}\\ {} {b^{-1}}\hspace{1em}& \text{if}\hspace{2.5pt}|b|\gt 1\text{.}\end{array}\right.\end{aligned}\end{array}\]
The desired formulas for the autocovariance function follow from (16) and (18).
The causality conditions, which are obtained in Proposition 4, can be rewritten as follows:
\[\begin{aligned}{}(1-b)(1-c)& \gt 0,\hspace{2em}(1+b)(1+c)\gt 0,\\ {} (1-b)(1+c)& \gt 0,\hspace{2em}(1+b)(1-c)\gt 0.\end{aligned}\]
These conditions are equivalent to $|b|\lt 1$ and $|c|\lt 1$. The proof is similar to one in Proposition 10.  □
The values of the autocovariance function for $a=-0.1$, $b=0.5$, $c=0.2$, and ${\sigma _{\epsilon }^{2}}=0.72$ are shown in Table 1. The parameters are chosen in such a way that $\operatorname{var}{X_{i,j}}=1$.
Table 1.
The autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$ for $a=-0.1$, $b=0.5$, $c=0.2$, and ${\sigma _{\epsilon }^{2}}=0.72$
Values of ${\boldsymbol{\gamma }_{\boldsymbol{X}}}\boldsymbol{(}{\boldsymbol{h}_{\mathbf{1}}}\mathbf{,}{\boldsymbol{h}_{\mathbf{2}}}\boldsymbol{)}$
${\boldsymbol{h}_{\mathbf{2}}}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\boldsymbol{-}\mathbf{2}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\boldsymbol{-}\mathbf{1}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\mathbf{0}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\mathbf{1}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\mathbf{2}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\mathbf{3}$
3 0 0 0.125 0.1125 0.0225 –0.00225
2 0 0 0.25 0.15 0.0075 –0.003
1 0 0 0.5 0.15 –0.015 0.0015
0 0 0 1 0 0 0
–1 –0.015 0.15 0.5 0 0 0
–2 0.0075 0.15 0.25 0 0 0
–3 0.0225 0.1125 0.125 0 0 0

5.2 Nonidentifiability. Symmetry

Let X be a stationary solution to equation (1). Consider the relation between parameters a, b, c and ${\sigma _{\epsilon }^{2}}$ and the distribution of the random field X – more specifically, we study the multiplicity of this relation. The situation is complicated by the fact that the distribution of X is not determined uniquely by the parameters a, b, c and ${\sigma _{\epsilon }^{2}}$; it also depends on the distribution of the error field ϵ besides the scaling coefficient ${\sigma _{\epsilon }}$.
Definition 2.
Let us have a statistical model $(\Omega ,\mathcal{F},\hspace{0.2778em}{\mathsf{P}_{\theta ,\phi }},\hspace{0.2222em}(\theta ,\phi )\in \Psi )$, that is, a family of probability measures ${\mathsf{P}_{\theta ,\psi }}$ on a common measurable space $(\Omega ,\mathcal{F})$ is indexed by parameter $(\theta ,\phi )\in \Psi $. The parameter θ is called identifiable if ${\mathsf{P}_{{\theta _{1}},{\psi _{1}}}}={\mathsf{P}_{{\theta _{2}},{\psi _{2}}}}$ implies ${\theta _{1}}={\theta _{2}}$.
In the model considered, the parameter of interest θ comprises the coefficients and the error variance, $\theta =(a,b,c,{\sigma _{\epsilon }^{2}})\in \Theta $, ${\sigma _{\epsilon }}\gt 0$,
\[\begin{aligned}{}\Theta =& \{(a,b,c,v):(1-a-b-c)(1-a+b+c)(1+a-b+c)(1+a+b-c)\gt 0\\ {} & \hspace{2.5pt}\text{and}\hspace{2.5pt}v\gt 0\}.\end{aligned}\]
The nuisance “nonparametric” parameter ϕ describes the distribution of the normalized error field ${\sigma _{\epsilon }^{-1}}\epsilon $. In the Gaussian case, ${\sigma _{\epsilon }^{-1}}\epsilon $ is a collection of independent random variables all of standard normal distribution. In general case, ${\sigma _{\epsilon }^{-1}}\epsilon $ is a collection of uncorrelated zero-mean unit-variance random variables.
It is much simpler to consider the autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$ than the distribution of X, as ${\gamma _{X}}({h_{1}},{h_{2}})$ is uniquely determined by θ due to (12). In the Gaussian case, the distribution of zero-mean stationary field X is uniquely determined by its autocovariance function; and obviously the distribution uniquely determines the autocovariance function. In general case, if parameters ${\theta _{1}}$ and ${\theta _{2}}$ determine the same autocovariance function, then for any “distribution of a collection of uncorrelated zero-mean unit-variance variables ${\sigma _{\epsilon }^{-1}}\epsilon $” described by parameter ${\phi _{1}}$ there exists a corresponding distribution described by parameter ${\phi _{2}}$ such that ${P_{{\theta _{1}},{\phi _{1}}}}={P_{{\theta _{2}},{\phi _{2}}}}$. Thus, it is enough to find out which parameters θ yield the same autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$.
In this section we will abuse notation by denoting the random field X as ${X_{i,j}}$. Similarly, ${X_{-i,-j}}$ will be a random field ${X^{(2)}}$ defined so that ${X_{i,j}^{(2)}}={X_{-i,-j}}$.
By substitution of indices, (1) can be rewritten as
\[\begin{array}{l}\displaystyle {X_{-i,-j}}=-\frac{b}{c}{X_{1-i,-j}}-\frac{a}{c}{X_{-i,1-j}}+\frac{1}{c}{X_{1-i,1-j}}-\frac{1}{c}{\epsilon _{1-i,1-j}},\\ {} \displaystyle {X_{-i,j}}=\frac{1}{a}{X_{1-i,j}}-\frac{c}{a}{X_{-1,j-1}}-\frac{b}{a}{X_{1-i,j-1}}-\frac{1}{a}{\epsilon _{1-i,j}},\\ {} \displaystyle {X_{i,-j}}=-\frac{c}{b}{X_{i-1,-j}}+\frac{1}{b}{X_{i,1-j}}-\frac{a}{b}{X_{i-1,1-j}}-\frac{1}{b}{\epsilon _{i,1-j}}.\end{array}\]
Thus, random fields ${X_{-i,-j}}$, ${X_{-i,j}}$ and ${X_{i,-j}}$ satisfy equations of the structure of (1), but with different coefficients and different error term. The correspondence of parameters is presented in Table 2.
However, the random fields ${X_{i,j}}$ and ${X_{-i,-j}}$ have the same autocovariance function, ${\gamma _{{X_{i,j}}}}({h_{1}},{h_{2}})={\gamma _{{X_{-i,-j}}}}({h_{1}},{h_{2}})$. Hence, two sets of parameters ${\theta _{1}}=(a,b,c,{\sigma _{\epsilon }^{2}})$ and ${\theta _{2}}=(-b/c,-a/c,1/c,{c^{-2}}{\sigma _{\epsilon }^{2}})$ determine the same autocovariance function of the field X.
Table 2.
Correspondence of terms of the recurrence equation for fields ${X_{i,j}}$, ${X_{-i,-j}}$, ${X_{-i,j}}$ and ${X_{i,-j}}$
Line Field of Coefficients Error field
no., $\boldsymbol{m}$ interest variance values
1 ${X_{i,j}}$ a b c ${\sigma _{\epsilon }^{2}}$ ${\epsilon _{i,j}}$
2 ${X_{-i,-j}}$ $-b/c$ $-a/c$ $1/c$ ${c^{-2}}{\sigma _{\epsilon }^{2}}$ $-{c^{-1}}{\epsilon _{1-i,1-j}}$
3 ${X_{-i,j}}$ $1/a$ $-c/a$ $-b/a$ ${a^{-2}}{\sigma _{\epsilon }^{2}}$ $-{a^{-1}}{\epsilon _{1-i,j}}$
4 ${X_{i,j-1}}$ $-c/b$ $1/b$ $-a/b$ ${b^{-2}}{\sigma _{\epsilon }^{2}}$ $-{b^{-1}}{\epsilon _{i,1-j}}$
The autocovariance function of random fields ${X_{-i,j}}$ and ${X_{i,-j}}$ is flip-symmetric to the autocovariance function of the field ${X_{i,j}}$:
\[ {\gamma _{{X_{-i,j}}}}({h_{1}},{h_{2}})={\gamma _{{X_{i,-j}}}}({h_{1}},{h_{2}})={\gamma _{{X_{i,j}}}}(-{h_{1}},{h_{2}})={\gamma _{{X_{i,j}}}}({h_{1}},-{h_{2}}).\]
In the symmetric case $ab+c=0$, where the autocovariance function ${\gamma _{{X_{i,j}}}}({h_{1}},{h_{2}})$ is even in each of the arguments, the autocovariance functions of the random fields ${X_{i,j}}$, ${X_{-i,-j}}$, ${X_{-i,j}}$ and ${X_{i,-j}}$ are all equal. Hence, four sets of parameters presented in Table 3 determine the same autocovariance function of the random fields ${X_{i,j}}$. Hence, four sets of parameters ${\theta _{1}}=(a,b,-bc,{\sigma _{\epsilon }^{2}})$, ${\theta _{2}}=({a^{-1}}\hspace{-0.1667em},{b^{-1}}\hspace{-0.1667em},-{a^{-1}}{b^{-1}}\hspace{-0.1667em},{a^{-2}}{b^{-2}}{\sigma _{\epsilon }^{2}})$, ${\theta _{3}}=({a^{-1}}\hspace{-0.1667em},b,-{a^{-1}}b,{a^{-2}}{\sigma _{\epsilon }^{2}})$ and ${\theta _{4}}=(a,{b^{-1}}\hspace{-0.1667em},-a{b^{-1}}\hspace{-0.1667em},{b^{-2}}{\sigma _{\epsilon }^{2}})$ determine the same autocovariance function of the random field X.
Table 3.
Four sets of parameters $\theta =(a,b,c,{\sigma _{\epsilon }^{2}})$ that determine the same autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$
a b $c=-ab$ ${\sigma _{\epsilon }^{2}}$
${a^{-1}}$ ${b^{-1}}$ $-{a^{-1}}{b^{-1}}$ ${a^{-2}}{b^{-2}}{\sigma _{\epsilon }^{2}}$
${a^{-1}}$ b $-{a^{-1}}b$ ${a^{-2}}{\sigma _{\epsilon }^{2}}$
a ${b^{-1}}$ $-a{b^{-1}}$ ${b^{-2}}{\sigma _{\epsilon }^{2}}$
Now we prove that there no other spurious cases where different sets of parameters determine the same autocovariance function.
Lemma 4.
Let $D\gt 0$, with D defined in Lemma 1. Then of four sets of parameters listed in Table 2, exactly one set satisfy the causality conditions specified in Proposition 4.
The proof involves the search of 7 cases of different signs of ${f_{1}},\dots ,{f_{4}}$ defined in (31), and evaluating inequalities. For brevity, we do not present the full proof; however, in Table 4 we list for what signs of ${f_{1}},\dots ,{f_{4}}$ which parameterization is causal.
Table 4.
Relation between the signs of ${f_{1}},\dots ,{f_{4}}$ and line in Table 2 where parameters satisfy the causality condition
Signs of ${\boldsymbol{f}_{\mathbf{1}}}\mathbf{,}\dots \mathbf{,}{\boldsymbol{f}_{\mathbf{4}}}$ The Which Outcome
signs parameter-
imply ization
that is causal, $\boldsymbol{m}$
${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$, ${f_{4}}\gt 0$ 1 ${X_{i,j}}$ is causal w.r.t. ${\epsilon _{i,j}}$
${f_{1}}\lt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$, ${f_{4}}\lt 0$ $c\gt 0$ 2 ${X_{-i,-j}}$ is causal
${f_{1}}\gt 0$, ${f_{2}}\lt 0$, ${f_{3}}\lt 0$, ${f_{4}}\gt 0$ $c\lt 0$ w.r.t. ${\epsilon _{1-i,1-j}}$
${f_{1}}\lt 0$, ${f_{2}}\lt 0$, ${f_{3}}\gt 0$, ${f_{4}}\gt 0$ $a\gt 0$ 3 ${X_{-i,j}}$ is causal
${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\lt 0$, ${f_{4}}\lt 0$ $a\lt 0$ w.r.t. ${\epsilon _{1-i,j}}$
${f_{1}}\lt 0$, ${f_{2}}\gt 0$, ${f_{3}}\lt 0$, ${f_{4}}\gt 0$ $b\gt 0$ 4 ${X_{i,-j}}$ is causal
${f_{1}}\gt 0$, ${f_{2}}\lt 0$, ${f_{3}}\gt 0$, ${f_{4}}\lt 0$ $b\lt 0$ w.r.t. ${\epsilon _{i,1-j}}$
Lemma 5.
In the causal case, that is, under causality conditions stated in Proposition 4, the autocovariance function uniquely determines coefficients a, b and c and error variance ${\sigma _{\epsilon }^{2}}$.
Proof.
Using the Yule–Walker equations, (20) and (33), we can express the parameters a, b, c and ${\sigma _{\epsilon }^{2}}$ in terms of ${\gamma _{X}}(0,0)$, ${\gamma _{X}}(1,0)$, ${\gamma _{X}}(0,1)$ and ${\gamma _{X}}(1,1)$. The explicit formulas are
\[\begin{array}{l}\displaystyle a=\frac{{\gamma _{X}}(1,0){\gamma _{X}}(0,0)-{\gamma _{X}}(0,1){\gamma _{X}}(1,1)}{{\gamma _{X}}{(0,0)^{2}}-{\gamma _{X}}{(0,1)^{2}}},\\ {} \displaystyle b=\frac{{\gamma _{X}}(0,1){\gamma _{X}}(0,0)-{\gamma _{X}}(1,0){\gamma _{X}}(1,1)}{{\gamma _{X}}{(0,0)^{2}}-{\gamma _{X}}{(1,0)^{2}}},\\ {} \displaystyle c=\frac{{\gamma _{X}}(1,1)-a{\gamma _{X}}(0,1)-b{\gamma _{X}}(1,0)}{{\gamma _{X}}(0,0)},\\ {} \displaystyle {\sigma _{\epsilon }^{2}}={\gamma _{X}}(0,0)\sqrt{D},\end{array}\]
where D is defined in Lemma 1.  □
The following proposition provides a necessary condition for two sets of parameters to determine the same autocovariance function of the field X.
Proposition 12.
Let stationary fields ${X^{(1)}}$ and ${X^{(2)}}$ satisfy (1)-like equation with parameters ${\theta _{1}}=({a_{1}},{b_{1}},{c_{1}},{\sigma _{\epsilon ,1}^{2}})$ and ${\theta _{2}}=({a_{2}},{b_{2}},{c_{2}},{\sigma _{\epsilon ,2}^{2}})$:
\[\begin{array}{l}\displaystyle {X_{i,j}^{(1)}}={a_{1}}{X_{i-1,j}^{(1)}}+{b_{1}}{X_{i,j-1}^{(1)}}+{c_{1}}{X_{i-1,j-1}^{(1)}}+{\epsilon _{i,j}^{(1)}},\\ {} \displaystyle {X_{i,j}^{(2)}}={a_{2}}{X_{i-1,j}^{(2)}}+{b_{2}}{X_{i,j-1}^{(2)}}+{c_{2}}{X_{i-1,j-1}^{(2)}}+{\epsilon _{i,j}^{(2)}},\\ {} \displaystyle \operatorname{var}{\epsilon _{i,j}^{(1)}}={\sigma _{\epsilon ,1}^{2}},\hspace{2em}\operatorname{var}{\epsilon _{i,j}^{(2)}}={\sigma _{\epsilon ,2}^{2}}.\end{array}\]
If random fields ${X^{(1)}}$ and ${X^{(2)}}$ have the same autocovariance function, then the parameters ${\theta _{1}}$ and ${\theta _{2}}$ relate to each other as quadruples of parameters in Table 2.
Proof.
Denote ${T_{m}}:{\Theta _{m}}\to \Theta $ the operator that transforms the quadruple of parameters in the fist line of Table 2 into one in the mth line, $m=1,2,3,4$. Thus,
\[\begin{array}{l}\displaystyle {T_{2}}(a,b,c,{\sigma _{\epsilon }^{2}})=\left(-\frac{b}{c},-\frac{a}{c},\frac{1}{c},{c^{-2}}{\sigma _{\epsilon }^{2}}\right),\\ {} \displaystyle {T_{3}}(a,b,c,{\sigma _{\epsilon }^{2}})=\left(\frac{1}{a},-\frac{c}{a},-\frac{b}{a},{a^{-2}}{\sigma _{\epsilon }^{2}}\right),\\ {} \displaystyle {T_{4}}(a,b,c,{\sigma _{\epsilon }^{2}})=\left(-\frac{c}{b},\frac{1}{b},-\frac{a}{b},{b^{-2}}{\sigma _{\epsilon }^{2}}\right),\end{array}\]
and ${T_{1}}$ is the identity operator, ${T_{1}}(\theta )=\theta $. Here ${\Theta _{m}}\subset \Theta $ is the domain where the operator ${T_{m}}$ is well defined, e.g., ${\Theta _{2}}=\{(a,b,c,{\sigma _{\epsilon }^{2}})\in \Theta :c\ne 0\}$. We have to prove that ${\theta _{1}}={T_{m}}{\theta _{0}}$ and ${\theta _{2}}={T_{n}}{\theta _{0}}$ for some $m,n=1,\dots ,4$ and ${\theta _{0}}\in \Theta $.
According to Lemma 4, one of parameterization ${T_{1}}{\theta _{1}},\dots ,{T_{4}}{\theta _{1}}$ satisfies the conditions for causality; let it be ${T_{m}}{\theta _{1}}={\theta _{3}}$. Denote also by ${X^{(3)}}$ and ${\epsilon ^{(3)}}$ the random field defined by such parameterization, and the respective white noise; what this means is listed in Table 2. For example, if $m=1$, that ${X^{(3)}}={X^{(1)}}$, and if $m=2$, then ${X_{i,j}^{(3)}}={X_{-i,-j}^{(1)}}$.
The random fields ${X^{(1)}}$ and ${X^{(3)}}$ are either equal or flip or turn symmetric to each other, ${X_{i,j}^{(3)}}={X_{\pm i,\pm j}^{(1)}}$. Hence, their autocovariance functions are either equal or one-variable symmetric to each other: either
\[ {\gamma _{{X^{(3)}}}}({h_{1}},{h_{2}})={\gamma _{{X^{(1)}}}}({h_{1}},{h_{2}})\hspace{1em}\text{for all}\hspace{2.5pt}{h_{1}}\hspace{2.5pt}\text{and}\hspace{2.5pt}{h_{2}},\]
or
\[ {\gamma _{{X^{(3)}}}}({h_{1}},{h_{2}})={\gamma _{{X^{(1)}}}}(-{h_{1}},{h_{2}})={\gamma _{{X^{(1)}}}}({h_{1}},-{h_{2}})\hspace{1em}\text{for all}\hspace{2.5pt}{h_{1}}\hspace{2.5pt}\text{and}\hspace{2.5pt}{h_{2}}.\]
The random field ${X^{(3)}}$ is causal w.r.t. white noise ${\epsilon ^{(3)}}$. Hence, due to Proposition 7,
(38)
\[ {\gamma _{{X^{(3)}}}}({h_{1}},{h_{2}})=\frac{{\gamma _{{X^{(3)}}}}({h_{1}},0){\gamma _{{X^{(3)}}}}(0,{h_{2}})}{{\gamma _{{X^{(3)}}}}(0,0)}\hspace{1em}\text{if}\hspace{2.5pt}{h_{1}}{h_{2}}\lt 0\text{.}\]
We do the same with the field ${X^{(2)}}$. There is n and a field ${X^{(4)}}$ that satisfy a (1)-like equation with parameters ${\theta _{4}}={T_{n}}{\theta _{2}}$; these parameters satisfy the conditions for causality, and either
\[\begin{array}{l}\displaystyle {\gamma _{{X^{(4)}}}}({h_{1}},{h_{2}})={\gamma _{{X^{(2)}}}}({h_{1}},{h_{2}})\hspace{1em}\text{for all}\hspace{2.5pt}{h_{1}}\hspace{2.5pt}\text{and}\hspace{2.5pt}{h_{2}}\text{, or}\\ {} \displaystyle {\gamma _{{X^{(4)}}}}({h_{1}},{h_{2}})={\gamma _{{X^{(2)}}}}(-{h_{1}},{h_{2}})={\gamma _{{X^{(2)}}}}({h_{1}},-{h_{2}})\hspace{1em}\text{for all}\hspace{2.5pt}{h_{1}}\hspace{2.5pt}\text{and}\hspace{2.5pt}{h_{2}},\end{array}\]
and also
(39)
\[ {\gamma _{{X^{(4)}}}}({h_{1}},{h_{2}})=\frac{{\gamma _{{X^{(4)}}}}({h_{1}},0){\gamma _{{X^{(4)}}}}(0,{h_{2}})}{{\gamma _{{X^{(4)}}}}(0,0)}\hspace{1em}\text{if}\hspace{2.5pt}{h_{1}}{h_{2}}\lt 0\text{.}\]
Similar relations hold for autocovariance functions of ${X^{(3)}}$ and ${X^{(4)}}$: either
\[\begin{array}{l}\displaystyle {\gamma _{{X^{(4)}}}}({h_{1}},{h_{2}})={\gamma _{{X^{(3)}}}}({h_{1}},{h_{2}})\hspace{1em}\text{for all}\hspace{2.5pt}{h_{1}}\hspace{2.5pt}\text{and}\hspace{2.5pt}{h_{2}}\text{, or}\\ {} \displaystyle {\gamma _{{X^{(4)}}}}({h_{1}},{h_{2}})={\gamma _{{X^{(3)}}}}(-{h_{1}},{h_{2}})={\gamma _{{X^{(3)}}}}({h_{1}},-{h_{2}})\hspace{1em}\text{for all}\hspace{2.5pt}{h_{1}}\hspace{2.5pt}\text{and}\hspace{2.5pt}{h_{2}},\end{array}\]
This implies that ${\gamma _{{X^{(4)}}}}({h_{1}},{h_{2}})={\gamma _{{X^{(3)}}}}({h_{1}},{h_{2}})$ if ${h_{1}}=0$ or ${h_{2}}=0$, and because of (38) and (39),
\[ {\gamma _{{X^{(4)}}}}({h_{1}},{h_{2}})={\gamma _{{X^{(3)}}}}({h_{1}},{h_{2}})\hspace{1em}\text{if}\hspace{2.5pt}{h_{1}}{h_{2}}\lt 0\text{.}\]
All above implies that random fields ${X^{(3)}}$ and ${X^{(4)}}$ have the same autocovariance function. These fields are causal w.r.t. respective white noises. According to Lemma 5, their parameters are equal, ${\theta _{3}}={\theta _{4}}$.
The operators ${T_{1}},\dots ,{T_{4}}$ are self-inverse; ${T_{k}^{2}}$ is the identity operator. Thus, ${T_{m}}{\theta _{1}}={\theta _{3}}={\theta _{4}}={T_{n}}{\theta _{2}}$ implies ${\theta _{1}}={T_{n}}{\theta _{3}}$ and ${\theta _{2}}={T_{m}}{\theta _{4}}$. Thus, parameters ${\theta _{1}}$ and ${\theta _{2}}$ relate to each other as the mth and nth rows of Table 2.  □
The following corollary combines necessary and sufficient conditions for parameters to determine the same autocovariance function; it also takes into account where the parameters make sense.
Corollary 2.
The parameter set Θ splits into two-element, four-element and one-element classes of parameters that define the same autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$ as follows:
  • 1. In generic case $-ab\ne c\ne 0$, the parameter $(a,b,c,{\sigma _{\epsilon }^{2}})$ belongs to a two-element class. It determines the same autocovariance function as the parameter $(-b/c,-a/c,1/c,{c^{-2}}{\sigma _{\epsilon }^{2}})$.
  • 2. In symmetric case $-ab=c\ne 0$, the parameter $(a,b,-bc,{\sigma _{\epsilon }^{2}})$ belongs to a four-element class. All four sets of parameters listed in Table 3 determine the same autocovariance function.
  • 3. In symmetric case $a\ne 0$, $b=c=0$, the parameter $(a,0,0,{\sigma _{\epsilon }^{2}})$ belongs to a two-element class. It determines the same autocovariance function as the parameter $({a^{-1}},0,0,{a^{-2}}{\sigma _{\epsilon }^{2}})$.
    Similarly, if $b\ne 0$, then the parameter $(0,b,0,{\sigma _{\epsilon }^{2}})$ belongs to two-element class and determines the same autocovariance function as the parameter $(0,{b^{-1}},0,{b^{-2}}{\sigma _{\epsilon }^{2}})$.
  • 4. For $-ab\ne c=0$, the parameter $(a,b,0,{\sigma _{\epsilon }^{2}})$ makes a class of its own. Likewise, for $a=b=c=0$, the parameter $(0,0,0,{\sigma _{\epsilon }^{2}})$ makes a class of its own. These are exceptional cases where the autocovariance function uniquely determines the parameters.
Next proposition shows how the asymmetric case is split into two major subcases.
Proposition 13.
Let X be a stationary field and ϵ be a collection of zero-mean equal-variance random variables that satisfy (1).
  • 1. If $ab+c=0$, then
    \[ {\gamma _{X}}({h_{1}},{h_{2}})=\frac{{\gamma _{X}}({h_{1}},0){\gamma _{X}}(0,{h_{2}})}{{\gamma _{X}}(0,0)}\hspace{1em}\textit{for all}\hspace{2.5pt}{h_{1}}\hspace{2.5pt}\textit{and}\hspace{2.5pt}{h_{2}}.\]
  • 2. If $ab+c\ne 0$ and $1+{c^{2}}\gt {a^{2}}+{b^{2}}$, then
    \[\begin{array}{l}\displaystyle {\gamma _{X}}(1,-1)=\frac{{\gamma _{X}}(1,0){\gamma _{X}}(0,1)}{{\gamma _{X}}(0,0)}\ne {\gamma _{X}}(1,1),\\ {} \displaystyle {\gamma _{X}}({h_{1}},{h_{2}})=\frac{{\gamma _{X}}({h_{1}},0){\gamma _{X}}(0,{h_{2}})}{{\gamma _{X}}(0,0)}\hspace{1em}\textit{for all}\hspace{2.5pt}{h_{1}}\hspace{2.5pt}\textit{and}\hspace{2.5pt}{h_{2}}\hspace{2.5pt}\textit{such that}\hspace{2.5pt}{h_{1}}{h_{2}}\le 0.\end{array}\]
  • 3. If $ab+c\ne 0$ and $1+{c^{2}}\lt {a^{2}}+{b^{2}}$, then
    \[\begin{array}{l}\displaystyle {\gamma _{X}}(1,-1)\ne \frac{{\gamma _{X}}(1,0){\gamma _{X}}(0,1)}{{\gamma _{X}}(0,0)}={\gamma _{X}}(1,1),\\ {} \displaystyle {\gamma _{X}}({h_{1}},{h_{2}})=\frac{{\gamma _{X}}({h_{1}},0){\gamma _{X}}(0,{h_{2}})}{{\gamma _{X}}(0,0)}\hspace{1em}\textit{for all}\hspace{2.5pt}{h_{1}}\hspace{2.5pt}\textit{and}\hspace{2.5pt}{h_{2}}\hspace{2.5pt}\textit{such that}\hspace{2.5pt}{h_{1}}{h_{2}}\ge 0.\end{array}\]
Proof.
The symmetric case $ab+c=0$ has been studied in Section 5.1.1. The desired equality follows from Proposition 10.
Let us borrow some notation from the proof of Proposition 12, with ${X^{(1)}}=X$. Let ${\theta _{3}}=({a_{3}},{b_{3}},{c_{3}},{\sigma _{\epsilon ,3}^{2}})={T_{m}}(a,b,c,{\sigma _{\epsilon }^{2}})$ be a set of parameters obtained from Lemma 4 that satisfy the causality conditions and let ${X^{(3)}}$ be a stationary random field defined by this parameterization; ${X^{(3)}}$ is causal w.r.t. the respective white noise. The sign of
\[ 1+{c^{2}}-{a^{2}}-{b^{2}}=\frac{{f_{1}}{f_{4}}+{f_{2}}{f_{3}}}{2}\]
determines the relation between the autocovariance functions of fields X and ${X^{(3)}}$. If $1+{c^{2}}\gt {a^{2}}+{b^{2}}$, then either ${X^{(3)}}=X$ or ${X_{i,j}^{(3)}}={X_{-i,-j}}$; in either case the fields X and ${X^{(3)}}$ have the same autocovariance function. Otherwise, if $1+{c^{2}}\lt {a^{2}}+{b^{2}}$, then either ${X_{i,j}^{(3)}}={X_{-i,j}}$ or ${X_{i,j}^{(3)}}={X_{i,-j}}$; in both cases the autocovariance functions are one-variable symmetric, ${\gamma _{{X^{(3)}}}}({h_{1}},{h_{2}})={\gamma _{X}}({h_{1}},-{h_{2}})$.
As ${X^{(3)}}$ is causal, due to Proposition 7
\[ {\gamma _{{X^{(3)}}}}({h_{1}},{h_{2}})=\frac{{\gamma _{{X^{(3)}}}}({h_{1}},0){\gamma _{{X^{(3)}}}}(0,{h_{2}})}{{\gamma _{{X^{(3)}}}}(0,0)}\hspace{1em}\text{if}\hspace{2.5pt}{h_{1}}{h_{2}}\le 0.\]
Hence, and from the relation between ${\gamma _{X}}$ and ${\gamma _{{X^{(3)}}}}$ all the desired equalities follow.
Now prove the inequalities. Assume that the would-be inequality is actually an equality, that is ${\gamma _{X}}(1,1)={\gamma _{X}}(1,0){\gamma _{X}}(0,1)/{\gamma _{X}}(0,0)$ in case $1+{c^{2}}\gt {a^{2}}+{b^{2}}$ or ${\gamma _{X}}(1,-1)={\gamma _{X}}(1,0){\gamma _{X}}(0,1)/{\gamma _{X}}(0,0)$ in case $1+{c^{2}}\lt {a^{2}}+{b^{2}}$. Then the autocovariance function of the field ${X^{(3)}}$ satisfies
\[ {\gamma _{{X^{(3)}}}}(1,1)=\frac{{\gamma _{{X^{(3)}}}}(1,0){\gamma _{{X^{(3)}}}}(0,1)}{{\gamma _{{X^{(3)}}}}(0,0)}.\]
The field ${X^{(3)}}$ is causal, and formulas in Lemma 5 yield
\[ {a_{3}}=\frac{{\gamma _{{X^{(3)}}}}(1,0)}{{\gamma _{{X^{(3)}}}}(0,0)},\hspace{1em}{b_{3}}=\frac{{\gamma _{{X^{(3)}}}}(0,1)}{{\gamma _{{X^{(3)}}}}(0,0)},\hspace{1em}{c_{3}}=-\frac{{\gamma _{{X^{(3)}}}}(1,0){\gamma _{{X^{(3)}}}}(0,1)}{{\gamma _{{X^{(3)}}}}{(0,0)^{2}}}=-{a_{3}}{b_{3}}.\]
Examining 4 cases, we can verify a similar relation for the original parameterization, $c=-ab$.  □

6 Pure nondeterminism

6.1 Sufficient condition

We borrow the definition of a pure nondeterministic random field from [16]. For the field on a planar lattice, the definition rewrites as follows, in term of subspaces of the Hilbert space of finite-variance random variables ${L^{2}}(\Omega ,\mathcal{F},\mathsf{P})$.
Definition 3.
Let X be a field of random variables of finite variances. X is called purely nondeterministic if and only if
(40)
\[\begin{aligned}{}\underset{i,j}{\overline{\text{span}}}{X_{i,j}}=\underset{i,j}{\overline{\text{span}}}\Big(& \overline{\text{span}}\{{X_{r,s}}:r\le i\}\cap {(\overline{\text{span}}\{{X_{r,s}}:r\lt i\})^{\perp }}\\ {} & \hspace{1em}\cap \overline{\text{span}}\{{X_{r,s}}:s\le j\}\cap {(\overline{\text{span}}\{{X_{r,s}}:s\lt j\})^{\perp }}\Big).\end{aligned}\]
Here ${\overline{\text{span}}_{i,j}}{X_{i,j}}$ is the smallest (closed) subspace of ${L^{2}}(\Omega ,\mathcal{F},\mathsf{P})$ that contains all observations of the field X, while $\overline{\text{span}}\{{X_{r,s}}:r\le i\}$ is a similar minimal subspace which contains all observations ${X_{r,s}}$ of the field X with the first index $r\le i$. The outer $\overline{\text{span}}$ in the right-hand side of (40) is the smallest subspace that contains all subspaces used as an argument.
Pure nondeterminism is a sufficient condition for a stationary field to be representable in form (29), where ϵ is some collection of uncorrelated variables, and coefficients ${\psi _{k,l}}$ satisfy ${\textstyle\sum _{k=0}^{\infty }}{\textstyle\sum _{l=0}^{\infty }}{\psi _{k,l}^{2}}\lt \infty $. This follows from [16, Theorem 2.1], while the necessity is easy to verify.
Proposition 14.
Let X be a stationary field and ϵ be a collection of uncorrelated zero-mean unit-variance random variables that satisfy (1). If $1+{c^{2}}\gt {a^{2}}+{b^{2}}$, then X is purely nondeterministic.
Proof.
We use notation D from Lemma 1 and notation ${f_{1}},\dots ,{f_{4}}$ from (31).
The existence of a stationary solution to (31) imply that $D\gt 0$. As shown in the proof of Proposition 13, inequalities $D\gt 0$ and $1+{c^{2}}\gt {a^{2}}+{b^{2}}$ hold true in three cases, ${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$, ${f_{4}}\gt 0$ or ${f_{1}}\lt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$, ${f_{4}}\lt 0$ or ${f_{1}}\gt 0$, ${f_{2}}\lt 0$, ${f_{3}}\lt 0$, ${f_{3}}\lt 0$.
If ${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$, ${f_{4}}\gt 0$, then the field X is causal with respect to white noise ϵ. Equations (1) and (29) imply that
\[ \overline{\text{span}}\{{X_{r,s}}:r\le i\}=\overline{\text{span}}\{{\epsilon _{r,s}}:r\le i\},\hspace{1em}\hspace{0.2778em}\overline{\text{span}}\{{X_{r,s}}:s\le j\}=\overline{\text{span}}\{{\epsilon _{r,s}}:s\le j\},\]
whence
\[\begin{array}{c}\displaystyle \overline{\text{span}}\{{X_{r,s}}:r\le i\}\cap {(\overline{\text{span}}\{{X_{r,s}}:r\lt i\})^{\perp }}\\ {} \displaystyle \cap \overline{\text{span}}\{{X_{r,s}}:s\le j\}\cap {(\overline{\text{span}}\{{X_{r,s}}:s\lt j\})^{\perp }}=\operatorname{span}{\epsilon _{i,j}}.\end{array}\]
Thus,
\[\begin{array}{c}\displaystyle \underset{i,j}{\overline{\text{span}}}\Big(\overline{\text{span}}\{{X_{r,s}}:r\le i\}\cap {(\overline{\text{span}}\{{X_{r,s}}:r\lt i\})^{\perp }}\\ {} \displaystyle \cap \overline{\text{span}}\{{X_{r,s}}:s\le j\}\cap {(\overline{\text{span}}\{{X_{r,s}}:s\lt j\})^{\perp }}\Big)=\underset{i,j}{\overline{\text{span}}}{\epsilon _{i,j}}=\underset{i,j}{\overline{\text{span}}}{X_{i,j}},\end{array}\]
and (40) holds true.
Now consider two other cases, where $-{f_{1}}$, ${f_{2}}$, ${f_{3}}$ and $-{f_{4}}$ are nonzero and of the same sign. In these cases $c\ne 0$, and
\[ {\tilde{\epsilon }_{i,j}}={X_{i,j}}+\frac{b}{c}{X_{i-1,j}}+\frac{a}{c}{X_{i,j-1}}-\frac{1}{c}{X_{i-1,j-1}}\]
is also a collection of uncorrelated zero-mean equal-variance random variables (to verify this, one can compute the spectral density of $\tilde{\epsilon }$). Causality condition in Proposition 4 can be easily verified; the stationary field X is causal w.r.t. white noise $\tilde{\epsilon }$. The field X is causal in these cases likewise.  □

6.2 Counterexample to Tjøstheim

It would be tempting to use the sufficient condition for pure nondeterminism from [16, Theorem 3.1]; however, as noted in [18], that condition is not correct. A counterexample is constructed below.
Let $0\lt |\theta |\lt 1$. Let X be a stationary field that satisfies the equation
\[ {X_{i-1,j}}=\theta {X_{i,j-1}}+{\epsilon _{i,j}},\]
where ϵ is a collection of independent random variables with standard normal distribution, ${\epsilon _{i,j}}\sim \mathcal{N}(0,1)$.
The spectral density of the field X is
\[ {f_{X}}({\nu _{1}},{\nu _{2}})=\frac{1}{|{e^{2\pi \mathrm{i}{\nu _{1}}}}-\theta {e^{2\pi \mathrm{i}{\nu _{2}}}}{|^{2}}}\hspace{0.1667em}.\]
The denominator $|{e^{2\pi \mathrm{i}{\nu _{1}}}}-\theta {e^{2\pi \mathrm{i}{\nu _{2}}}}{|^{2}}$ attains only positive values and is a continuous function. Thus, the field X satisfies the sufficient condition stated in [16, Theorem 3.1].
In the field X the diagonals $\{{X_{i,j}},\hspace{0.2778em}i+j=n\}$ are independent for different n. The random variables on each diagonal are jointly distributed as values of a centered Gaussian AR(1) process with the coefficient θ.
Let
\[ {\tilde{\epsilon }_{i,j}}={X_{i,j-1}}-\theta {X_{i-1,j}}.\]
Then $\tilde{\epsilon }$ is a collection of uncorrelated zero-mean unit-variance variables. (Since $\tilde{\epsilon }$ is a Gaussian field, it is a collection independent variables with distribution $\mathcal{N}(0,1)$.)
The field X can be represented as
\[ {X_{i,j}}={\sum \limits_{k=0}^{\infty }}{\theta ^{k}}{\epsilon _{i+k+1,j-k}}={\sum \limits_{k=0}^{\infty }}{\theta ^{k}}{\tilde{\epsilon }_{i-k,j+k+1}}.\]
Hence,
\[ {\tilde{\epsilon }_{i,j}}=-\theta {\epsilon _{i,j}}+(1-{\theta ^{2}}){\sum \limits_{k=0}^{\infty }}{\epsilon _{i+k+1,j-k-1}}.\]
These representations imply that
\[ \overline{\text{span}}\{{X_{r,s}}:r\le i\}=\overline{\text{span}}\{{\tilde{\epsilon }_{r,s}}:r\le i\},\hspace{1em}\overline{\text{span}}\{{X_{r,s}}:s\le j\}=\overline{\text{span}}\{{\epsilon _{r,s}}:s\le j\}.\]
Hence,
\[\begin{array}{c}\displaystyle \overline{\text{span}}\{{X_{r,s}}:r\le i\}\cap {(\overline{\text{span}}\{{X_{r,s}}:r\lt i\})^{\perp }}\\ {} \displaystyle \cap \overline{\text{span}}\{{X_{r,s}}:s\le j\}\cap {(\overline{\text{span}}\{{X_{r,s}}:s\lt j\})^{\perp }}=\underset{s}{\overline{\text{span}}}{\tilde{\epsilon }_{i,s}}\cap \underset{r}{\overline{\text{span}}}{\epsilon _{r,j}}.\end{array}\]
Now prove that ${\overline{\text{span}}_{s}}{\tilde{\epsilon }_{i,s}}\cap {\overline{\text{span}}_{r}}{\epsilon _{r,j}}$ is a trivial subspace. Let $\zeta \in {\overline{\text{span}}_{s}}{\tilde{\epsilon }_{i,s}}\cap {\overline{\text{span}}_{r}}{\epsilon _{r,j}}$,
\[ \zeta ={\sum \limits_{s=-\infty }^{\infty }}{k_{s}}{\tilde{\epsilon }_{i,s}}={\sum \limits_{r=-\infty }^{\infty }}{c_{r}}{\epsilon _{r,j}}.\]
Here the coefficients satisfy ${\textstyle\sum _{s=-\infty }^{\infty }}{k_{s}^{2}}={\textstyle\sum _{r=-\infty }^{\infty }}{c_{r}^{2}}\lt \infty $; the series converge in mean squares. The covariance between ${\tilde{\epsilon }_{i,s}}$ and ${\epsilon _{r,j}}$ is
\[ \operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,s}}{\epsilon _{r,j}}=\left\{\begin{array}{l@{\hskip10.0pt}l}-\theta \hspace{1em}& \text{if}\hspace{2.5pt}i=r\hspace{2.5pt}\text{and}\hspace{2.5pt}s=j\text{,}\\ {} {\theta ^{r-i-1}}(1-{\theta ^{2}})\hspace{1em}& \text{if}\hspace{2.5pt}r-i=s-j\gt 0\text{,}\\ {} 0\hspace{1em}& \text{otherwise.}\end{array}\right.\]
Anyway, $\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,s}}{\epsilon _{r,j}}=0$ if $r-i\ne s-j$, and $|\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,s}}{\epsilon _{r,j}}|\lt 1$ if $r-i=s-j$.
Compute $\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,j-n}}\zeta $ and $\operatorname{\mathsf{E}}\zeta {\epsilon _{i-n,j}}$ using two different expressions:
\[\begin{array}{l}\displaystyle {\sum \limits_{s=-\infty }^{\infty }}{k_{s}}\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,j-n}}{\tilde{\epsilon }_{i,s}}=\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,j-n}}\zeta ={\sum \limits_{r=-\infty }^{\infty }}{c_{r}}\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,j-n}}{\epsilon _{r,j}},\\ {} \displaystyle {\sum \limits_{s=-\infty }^{\infty }}{k_{s}}\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,s}}{\epsilon _{i-n,j}}=\operatorname{\mathsf{E}}\zeta {\epsilon _{i-n,j}}={\sum \limits_{r=-\infty }^{\infty }}{c_{r}}\operatorname{\mathsf{E}}{\epsilon _{r,j}}{\epsilon _{i-n,j}}.\end{array}\]
In the series, the only nonzero term might be where $s=j-n$ and $r=i-n$. Thus,
(41)
\[ {k_{j-n}}={c_{i-n}}\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,j-n}}{\epsilon _{i-n,j}},\hspace{2em}{k_{j-n}}\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,j-n}}{\epsilon _{i-n,j}}={c_{i-n}}.\]
As $|\operatorname{\mathsf{E}}{\tilde{\epsilon }_{i,j-n}}{\epsilon _{i-n,j}}|\lt 1$, (41) imply ${k_{j-n}}={c_{i-n}}=0$. As this holds true for all integer n, ${k_{s}}={c_{r}}=0$ for all s and r, and $\zeta =0$ almost surely. Thus,
\[\begin{array}{c}\displaystyle \overline{\text{span}}\{{X_{r,s}}:r\le i\}\cap {(\overline{\text{span}}\{{X_{r,s}}:r\lt i\})^{\perp }}\\ {} \displaystyle \cap \overline{\text{span}}\{{X_{r,s}}:s\le j\}\cap {(\overline{\text{span}}\{{X_{r,s}}:s\lt j\})^{\perp }}=\{0\},\end{array}\]
and the right-hand side of (40) is the trivial subspace. Thus, the random field X is not purely nondeterministic.

7 Conclusion

We considered AR(1) model on a plane. We found conditions (in terms of the regression coefficients) under which the autoregressive equation has a stationary solution X. As for the autocovariance function of the stationary solution X, we presented a simple formula for it at some points, and proved Yule–Walker equations. These allow to compute the autocovariance function recursively at all points.
We found conditions under which the stationary solution to the autoregressive equation satisfies the causality condition with respect to the underlying white noise. These conditions also appear to be sufficient conditions for stability of the deterministic problem of solving a recursive equations in a quadrant, with preset values on the border of the quadrant.
We described sets of parameters (the coefficients and the variance of the underlying white noise) where different parameters determine the same autocovariance function of the stationary solution.
We found sufficient conditions for the stationary solution X to be a pure nondeterministic random field. This condition is related to the causality condition with respect to some (nonfixed) white noise, which is called innovations; in particular, the innovations need not coincide with the white noise in the autoregressive equation.
The causality condition and pure-nondeterministic property seem to be too restrictive because the coordinate-wise order between coordinates of ϵ and X in the representation (29) is a partial order. More general representation with the lexical order, which is a total order, is suggested in [21, Section 6]:
\[ {X_{i,j}}={\sum \limits_{k=1}^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}{\psi _{k,l}}{X_{i-k,j-l}}+{\sum \limits_{l=0}^{\infty }}{\psi _{0,l}}{X_{i,j-l}}.\]
A further discussion can be found in [9, 17].
With exception for the identifiability topic, we did not consider the estimation at all.
Most results of this paper are necessary and sufficient conditions for some properties of the stationary field X. However, for pure nondeterminism of X and for stability of the deterministic equation, we obtained only sufficient conditions. Obtaining necessary conditions is an interesting open problem.

A Appendix: Auxiliary results

A.1 Existence of a random field with given spectral density

The next lemma is used in the proof of Proposition 1. It is a modification of the similar result stated for stochastic processes [15].
Lemma 6.
Let $f({\nu _{1}},{\nu _{2}})$ be an even integrable function ${\left[-\frac{1}{2},\frac{1}{2}\right]^{2}}\to [0,\infty ]$, that is,
\[\begin{array}{l}\displaystyle f({\nu _{1}},{\nu _{2}})=f(-{\nu _{1}},-{\nu _{2}})\ge 0\hspace{1em}\textit{for all}\hspace{2.5pt}\nu 1,{\nu _{2}}\in [-\frac{1}{2},\frac{1}{2}],\\ {} \displaystyle {\iint _{{[-1/2,\hspace{0.1667em}1/2]^{2}}}}f({\nu _{1}},{\nu _{2}})\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}\lt \infty .\end{array}\]
Then there exists a Gaussian stationary field $\{{X_{i,j}},\hspace{2.5pt}i,j\in \mathbb{Z}\}$ on some (specially constructed) probability space that has spectral density $f({\nu _{1}},{\nu _{2}})$.
Proof.
Let us construct a probability space with two independent Brownian fields on a plane $\{{W_{k}}(t,s),\hspace{2.5pt}t,s\in [0,1]\}$ (a zero-mean Gaussian field with covariance function $\operatorname{cov}({W_{k}}({t_{1}},{s_{1}}),\hspace{0.1667em}{W_{k}}({t_{2}},{s_{2}}))=\min ({t_{1}},{t_{2}})\min ({s_{1}},{s_{2}})$, $k=1,\hspace{0.1667em}2$).
Extend the function X by periodicity, $f({\nu _{1}},{\nu _{2}})=f({\nu _{1}}-1,{\nu _{2}})$ if $\frac{1}{2}\lt {\nu _{1}}\le 1$, and then $f({\nu _{1}},{\nu _{2}})=f({\nu _{1}},{\nu _{2}}-1)$ if $\frac{1}{2}\lt {\nu _{2}}\le 1$.
Denote
\[\begin{aligned}{}{X_{i,j}}& ={\iint _{{[0,1]^{2}}}}\cos (2\pi (i{\nu _{1}}+j{\nu _{2}}))\sqrt{f({\nu _{1}},{\nu _{2}})}\hspace{0.1667em}{d^{2}}{W_{1}}({\nu _{1}},{\nu _{2}})\\ {} & +{\iint _{{[0,1]^{2}}}}\sin (2\pi (i{\nu _{1}}+j{\nu _{2}}))\sqrt{f({\nu _{1}},{\nu _{2}})}\hspace{0.1667em}{d^{2}}{W_{2}}({\nu _{1}},{\nu _{2}}).\end{aligned}\]
The constructed field is zero-mean Gaussian as the integral of nonstochastic kernel with respect to the Gaussian field. The autocovariance function is as expected,
\[\begin{aligned}{}\operatorname{cov}({X_{i,j}},{X_{i+{h_{1}},j+{h_{2}}}})& ={\iint _{{[0,1]^{2}}}}\cos (2\pi ({h_{1}}{\nu _{1}}+{h_{2}}{\nu _{2}}))f({\nu _{1}},{\nu _{2}})\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}\\ {} & ={\iint _{{[-1/2,1/2]^{2}}}}\exp (2\pi \mathrm{i}({h_{1}}{\nu _{1}}+{h_{2}}{\nu _{2}}))f({\nu _{1}},{\nu _{2}})\hspace{0.1667em}d{\nu _{1}}\hspace{0.1667em}d{\nu _{2}}.\end{aligned}\]
 □

A.2 Summation

The next two propositions are used implicitly when we use double-series notation.
Proposition 15.
Let $\{{\xi _{i,j}},\hspace{0.2778em}i,j\in \mathbb{Z}\}$ be a collection of random variables bounded in mean squares, ${\sup _{i,j}}\operatorname{\mathsf{E}}{\xi _{i,j}^{2}}\lt \infty $. Let $\{{\xi _{i,j}},\hspace{0.2778em}i,j\in \mathbb{Z}\}$ be a collection of real numbers such that ${\textstyle\sum \textstyle\sum _{i,j=-\infty }^{\infty }}|{a_{i,j}}|\lt \infty $. Then the double series ${\textstyle\sum \textstyle\sum _{i,j=-\infty }^{\infty }}{a_{i,j}}{\xi _{i,j}}$ converges in mean squares and almost surely. The limit is the same for both types of convergence (up to equal-almost-surely equivalence); it also does not depend on the order the terms of the double series are added.
Remark 5.
  • 1. Under conditions of Proposition 15, the double series ${\textstyle\sum \textstyle\sum _{i,j=-\infty }^{\infty }}{a_{i,j}}{\xi _{i,j}}$ converges in mean and in probability, as well – to the same limit.
  • 2. The iterated series ${\textstyle\sum _{i=-\infty }^{\infty }}\left({\textstyle\sum _{j=-\infty }^{\infty }}{a_{i,j}}{\xi _{i,j}}\right)$ converges to the same limit.
  • 3. $\operatorname{\mathsf{P}}\left({\textstyle\sum \textstyle\sum _{i,j=-\infty }^{\infty }}|{a_{i,j}}{\xi _{i,j}}|\lt \infty \right)=1$, and whenever ${\textstyle\sum \textstyle\sum _{i,j=-\infty }^{\infty }}|{a_{i,j}}{\xi _{i,j}}|\lt \infty $, the double series ${\textstyle\sum \textstyle\sum _{i,j=-\infty }^{\infty }}{a_{i,j}}{\xi _{i,j}}$ converges, and its limit does not depend on the order its terms are added.
Proposition 16.
Let $\{{\epsilon _{i,j}},\hspace{0.2778em}i,j\in \mathbb{Z}\}$ be a collection of random uncorrelated variables with zero mean and same variance $\operatorname{var}{\epsilon _{i,j}}={\sigma _{\epsilon }^{2}}\lt \infty $. Let $\{{\xi _{i,j}},\hspace{0.2778em}i,j\in \mathbb{Z}\}$ be a collection of real numbers such that ${\textstyle\sum \textstyle\sum _{i,j=-\infty }^{\infty }}{a_{i,j}^{2}}\lt \infty $. Then the double series ${\textstyle\sum \textstyle\sum _{i,j=-\infty }^{\infty }}{a_{i,j}}{\xi _{i,j}}$ converges in mean squares. The limit does not depend on the order the terms of the double series are added: it is the same up to equal-almost-surely equivalence.
Remark 6.
  • 1. Under conditions of Proposition 16, the double series ${\textstyle\sum \textstyle\sum _{i,j=-\infty }^{\infty }}{a_{i,j}}{\xi _{i,j}}$ converges in mean and in probability, as well.
  • 2. The iterated series ${\textstyle\sum _{i=-\infty }^{\infty }}\left({\textstyle\sum _{j=-\infty }^{\infty }}{a_{i,j}}{\xi _{i,j}}\right)$ converges in mean squares to the same limit.
The next lemma is used in Proposition 4.
Lemma 7.
Let $f\in C({\mathbb{R}^{2}})$ be a biperiodic function of two variables with continuous mixed derivative:
\[\begin{array}{l}\displaystyle f(x,y)=f(x+1,y)=f(x,y+1)\hspace{1em}\textit{for all}\hspace{2.5pt}x,y\in \mathbb{R},\\ {} \displaystyle f\in C({\mathbb{R}^{2}}),\hspace{1em}\frac{\partial f}{\partial x}\in C({\mathbb{R}^{2}}),\hspace{1em}\frac{\partial f}{\partial y}\in C({\mathbb{R}^{2}}),\hspace{1em}\frac{{\partial ^{2}}f}{\partial x\hspace{0.1667em}\partial y}\in C({\mathbb{R}^{2}}).\end{array}\]
Then Fourier coefficients of the function f are summable:
\[ {\sum \limits_{k=-\infty }^{\infty }}{\sum \limits_{l=-\infty }^{\infty }}\left|{\int _{0}^{1}}{\int _{0}^{1}}\exp (2\pi \mathrm{i}\hspace{0.1667em}(kx+ly))\hspace{0.1667em}f(x,y)\hspace{0.1667em}dx\hspace{0.1667em}dy\right|\lt \infty .\]
Proof.
Denote the Fourier coefficients ${c_{k,l}}$:
\[ {c_{k,l}}={\int _{0}^{1}}{\int _{0}^{1}}\exp (2\pi \mathrm{i}\hspace{0.1667em}(kx+ly))\hspace{0.1667em}f(x,y)\hspace{0.1667em}dx\hspace{0.1667em}dy.\]
The coefficients ${c_{k,l}}$ are also equal
\[\begin{array}{l}\displaystyle {c_{k,l}}=\frac{\mathrm{i}}{2\pi k}{\int _{0}^{t}}{\int _{0}^{1}}\exp (2\pi \mathrm{i}\hspace{0.1667em}(kx+ly))\hspace{0.1667em}{f^{\prime }_{1}}(x,y)\hspace{0.1667em}dx\hspace{0.1667em}dy\hspace{1em}\text{if}\hspace{2.5pt}k\ne 0;\\ {} \displaystyle {c_{k,l}}=\frac{\mathrm{i}}{2\pi l}{\int _{0}^{t}}{\int _{0}^{1}}\exp (2\pi \mathrm{i}\hspace{0.1667em}(kx+ly))\hspace{0.1667em}{f^{\prime }_{2}}(x,y)\hspace{0.1667em}dx\hspace{0.1667em}dy\hspace{1em}\text{if}\hspace{2.5pt}l\ne 0;\\ {} \displaystyle {c_{k,l}}=\frac{-1}{4{\pi ^{2}}kl}{\int _{0}^{t}}{\int _{0}^{1}}\exp (2\pi \mathrm{i}\hspace{0.1667em}(kx+ly))\hspace{0.1667em}{f^{\prime\prime }_{12}}(x,y)\hspace{0.1667em}dx\hspace{0.1667em}dy\hspace{1em}\text{if}\hspace{2.5pt}k\ne 0\hspace{2.5pt}\text{and}\hspace{2.5pt}l\ne 0,\end{array}\]
where ${f^{\prime }_{1}}(x,y)$, ${f^{\prime }_{2}}(x,y)$, and ${f^{\prime\prime }_{12}}(x,y)$ are partial and mixed derivatives; here the periodicity of the function f is used. The coefficients ${c_{k,l}}$ allow bounding:
\[\begin{array}{l}\displaystyle {c_{k,0}}\le \frac{1}{2\pi |k|}\underset{x,y}{\max }\left|\frac{\partial f(x,y)}{\partial x}\right|\hspace{1em}\text{if}\hspace{2.5pt}k\ne 0;\\ {} \displaystyle {c_{0,l}}\le \frac{1}{2\pi |l|}\underset{x,y}{\max }\left|\frac{\partial f(x,y)}{\partial y}\right|\hspace{1em}\text{if}\hspace{2.5pt}l\ne 0;\\ {} \displaystyle {c_{k,l}}\le \frac{1}{4{\pi ^{2}}|kl|}\underset{x,y}{\max }\left|\frac{{\partial ^{2}}f(x,y)}{\partial y\partial x}\right|\hspace{1em}\text{if}\hspace{2.5pt}k\ne 0\hspace{2.5pt}\text{and}\hspace{2.5pt}l\ne 0,\end{array}\]
which implies the summability.  □

A.3 Integration

Lemma 8.
If $A\gt 0$, B and C are real numbers such that ${A^{2}}\gt {B^{2}}+{C^{2}}$, then
(42)
\[ {\int _{-1/2}^{1/2}}\frac{dt}{A+B\cos (2\pi t)+C\sin (2\pi t)}=\frac{1}{\sqrt{{A^{2}}-{B^{2}}-{C^{2}}}}.\]
If A and B are real numbers, $A\gt |B|$, and n is integer, then
(43)
\[ {\int _{-1/2}^{1/2}}\frac{\exp (2\pi \mathrm{i}nt)}{A+B\cos (2\pi t)}\hspace{0.1667em}dt=\frac{{\alpha ^{|n|}}}{\sqrt{{A^{2}}-{B^{2}}}},\]
where
\[ \alpha =\frac{-B}{A+\sqrt{{A^{2}}-{B^{2}}}}=\left\{\begin{array}{l@{\hskip10.0pt}l}0\hspace{1em}& \textit{if}\hspace{2.5pt}B=0,\\ {} (-A+\sqrt{{A^{2}}-{B^{2}}})/B\hspace{1em}& \textit{if}\hspace{2.5pt}B\ne 0.\end{array}\right.\]
In (43), if $\alpha =0$ and $n=0$, then ${\alpha ^{|n|}}=1$ by convention.
Proof.
That is easy to check that the antiderivative in (42) is
\[\begin{aligned}{}& \int \frac{dt}{A+B\cos (2\pi t)+C\sin (2\pi t)}\\ {} & \hspace{1em}=\frac{1}{\pi \sqrt{{A^{2}}-{B^{2}}-{C^{2}}}}\hspace{0.1667em}\arctan \hspace{-0.1667em}\left(\frac{(A-B)\tan (\pi t)+C}{\sqrt{{A^{2}}-{B^{2}}-{C^{2}}}}\right)+\mathrm{const},\end{aligned}\]
whence (42) follows. (Notice that $A-B\gt 0$.)
In (43), the function
\[ \frac{1}{A+B\cos (2\pi t)}=\frac{1+{\alpha ^{2}}}{A\hspace{0.1667em}(1-2\alpha \cos (2\pi t)+{\alpha ^{2}})}=\frac{1+{\alpha ^{2}}}{A}\hspace{0.1667em}\frac{1}{|1-\alpha \exp (2\pi \mathrm{i}t){|^{2}}}\]
is the spectral density of AR(1) stationary autoregressive process ${X_{k}}=\alpha {X_{k-1}}+{\epsilon _{k}}$ with the white noise variance $\operatorname{var}{\epsilon _{k}}=(1+{\alpha ^{2}})/A$. (Conditions imply that $|\alpha |\lt 1$.) The integral is the autocovariance function of the process. It equals
\[ {\int _{-1/2}^{1/2}}\frac{\exp (2\pi \mathrm{i}nt)}{A+\cos (2\pi t)}\hspace{0.1667em}dt=\operatorname{cov}({X_{k+n}},{X_{k}})=\frac{1+{\alpha ^{2}}}{A}\hspace{0.1667em}\frac{{\alpha ^{|n|}}}{1-{\alpha ^{2}}}=\frac{{\alpha ^{|n|}}}{\sqrt{{A^{2}}-{B^{2}}}}\hspace{0.1667em};\]
here equality
\[ \frac{1+{\alpha ^{2}}}{1-{\alpha ^{2}}}=\frac{A}{\sqrt{{A^{2}}-{B^{2}}}}\]
is used.  □
Lemma 9.
If A and B are complex numbers, $|A|\ne |B|$, and n is an integer number, then
\[ {\int _{-1/2}^{1/2}}\frac{{e^{2\pi \mathrm{i}n\nu }}}{A-B{e^{2\pi \mathrm{i}\nu }}}\hspace{0.1667em}d\nu =\left\{\begin{array}{l@{\hskip10.0pt}l}{A^{n-1}}{B^{-n}}\hspace{1em}& \textit{if}\hspace{2.5pt}n\le 0\hspace{2.5pt}\textit{and}\hspace{2.5pt}|A|\gt |B|\textit{,}\\ {} 0\hspace{1em}& \textit{if}\hspace{2.5pt}n\le 0\hspace{2.5pt}\textit{and}\hspace{2.5pt}|A|\lt |B|\textit{,}\\ {} 0\hspace{1em}& \textit{if}\hspace{2.5pt}n\ge 1\hspace{2.5pt}\textit{and}\hspace{2.5pt}|A|\gt |B|\textit{,}\\ {} -{A^{n-1}}{B^{-n}}\hspace{1em}& \textit{if}\hspace{2.5pt}n\ge 1\hspace{2.5pt}\textit{and}\hspace{2.5pt}|A|\lt |B|\textit{.}\end{array}\right.\]
Proof.
For $n=1$, the integral can be rewritten as the contour integral,
\[ {\int _{-1/2}^{1/2}}\frac{{e^{2\pi \mathrm{i}\nu }}}{A-B{e^{2\pi \mathrm{i}\nu }}}\hspace{0.1667em}d\nu =\frac{1}{2\pi \mathrm{i}}{\oint _{z=\exp (2\pi \mathrm{i}\nu )}}\frac{1}{A-Bz}\hspace{0.1667em}dz,\]
and the residue formula is applicable. Other cases can be reduced to the case $n=1$ recursively.  □
Lemma 10.
Let a, b and c be real numbers such that $1-a-b-c\gt 0$, $1-a+b+c\gt 0$, $1+a-b+c\gt 0$ and $1+a+b-c\gt 0$. Let ${n_{1}}$ and ${n_{2}}$ be nonnegative integers. Then
\[ {\int _{-1/2}^{1/2}}\frac{{(a+c{e^{2\pi \mathrm{i}\nu }})^{{n_{1}}}}{e^{-2\pi \mathrm{i}{n_{2}}\nu }}}{{(1-b{e^{2\pi \mathrm{i}\nu }})^{{n_{1}}+1}}}\hspace{0.1667em}d\nu ={\sum \limits_{k=0}^{\min ({n_{1}},{n_{2}})}}\left(\genfrac{}{}{0.0pt}{}{{n_{1}}}{k}\right)\left(\genfrac{}{}{0.0pt}{}{{n_{2}}}{k}\right){a^{{n_{1}}-k}}{b^{{n_{2}}-k}}{(ab+c)^{k}}.\]
Proof.
Conditions of the lemma imply that
\[\begin{array}{l}\displaystyle 1-b=\frac{(1-a-b-c)+(1+a-b+c)}{2}\gt 0,\\ {} \displaystyle 1+b=\frac{(1-a+b+c)+(1+a+b-c)}{2}\gt 0,\\ {} \displaystyle |b|\lt 1.\end{array}\]
Let $k\ge 0$ and n be integer numbers. By the binomial formula,
(44)
\[ \frac{1}{{(1-bz)^{k+1}}}={\sum \limits_{m=0}^{\infty }}\left(\genfrac{}{}{0.0pt}{}{-k-1}{m}\right){(-bz)^{m}}={\sum \limits_{m=0}^{\infty }}\left(\genfrac{}{}{0.0pt}{}{k+m}{m}\right){(bz)^{m}}\]
for all complex z such that $|bz|\lt 1$, in particular, for all z such that $|z|\le 1$. Here $\left(\genfrac{}{}{0.0pt}{}{-k-1}{m}\right)=\frac{(-k-1)\cdots (-k-m)}{1\cdots m}={(-1)^{m}}\left(\genfrac{}{}{0.0pt}{}{k+m}{m}\right)$ is a coefficient in a binomial series. The rational function ${z^{n}}{(1-bz)^{-k-1}}$ has a singularity at point $1/b$ and a potential singularity (if $n\lt 0$) at point 0. The point $1/b$ lies outside the unit circle, and point 0 lies inside the unit circle. The expansion (44) implies that
\[\begin{aligned}{}{\operatorname{Res}_{z=0}}\left(\frac{{z^{n}}}{{(1-bz)^{k+1}}}\right)& ={\operatorname{Res}_{z=0}}\left({\sum \limits_{m=0}^{\infty }}\left(\genfrac{}{}{0.0pt}{}{k+m}{m}\right){b^{m}}{z^{m+n}}\right)\\ {} & =\left\{\begin{array}{l@{\hskip10.0pt}l}0\hspace{1em}& \text{if}\hspace{2.5pt}n\ge 0,\\ {} \left(\genfrac{}{}{0.0pt}{}{k-n-1}{-n-1}\right){b^{-n-1}}\hspace{1em}& \text{if}\hspace{2.5pt}n\le -1\hspace{2.5pt}\text{.}\end{array}\right.\end{aligned}\]
By the residue formula,
\[\begin{aligned}{}\oint \frac{{z^{n}}}{{(1-bz)^{k+1}}}\hspace{0.1667em}dz& =2\pi \mathrm{i}{\operatorname{Res}_{z=0}}\left(\frac{{z^{n}}}{{(1-bz)^{k+1}}}\right)\\ {} & =\left\{\begin{array}{l@{\hskip10.0pt}l}0\hspace{1em}& \text{if}\hspace{2.5pt}n\ge 0,\\ {} 2\pi \mathrm{i}\left(\genfrac{}{}{0.0pt}{}{k-n-1}{-n-1}\right){b^{-n-1}}\hspace{1em}& \text{if}\hspace{2.5pt}n\le -1\text{,}\end{array}\right.\end{aligned}\]
where the contour integral is taken along the unit circle contour.
Recall that ${n_{2}}\ge 0$ is integer. Then
(45)
\[\begin{array}{cc}& \displaystyle {\int _{-1/2}^{1/2}}\frac{\exp (2\pi \mathrm{i}(k-{n_{2}}){\nu _{2}})}{{(1-b{e^{2\pi \mathrm{i}{\nu _{2}}}})^{1+k}}}\hspace{0.1667em}d{\nu _{2}}=\frac{1}{2\pi \mathrm{i}}\oint \frac{{z^{k-{n_{2}}-1}}}{{(1-bz)^{k+1}}}\hspace{0.1667em}dz\\ {} & \displaystyle =\left\{\begin{array}{l@{\hskip10.0pt}l}0\hspace{1em}& \text{if}\hspace{2.5pt}k\gt {n_{2}}\text{,}\\ {} \left(\genfrac{}{}{0.0pt}{}{{n_{2}}}{{n_{2}}-k}\right){b^{{n_{2}}-k}}\hspace{1em}& \text{if}\hspace{2.5pt}k\le {n_{2}}\end{array}\right.=\left\{\begin{array}{l@{\hskip10.0pt}l}0\hspace{1em}& \text{if}\hspace{2.5pt}k\gt {n_{2}}\text{,}\\ {} \left(\genfrac{}{}{0.0pt}{}{{n_{2}}}{k}\right){b^{{n_{2}}-k}}\hspace{1em}& \text{if}\hspace{2.5pt}k\le {n_{2}}\text{.}\end{array}\right.\end{array}\]
The constant ${n_{1}}\ge 0$ is also integer, and by the binomial formula
(46)
\[\begin{aligned}{}{(a+c{e^{2\pi \mathrm{i}\nu }})^{{n_{1}}}}& ={((ab+c){e^{2\pi \mathrm{i}\nu }}+a(1-b{e^{2\pi \mathrm{i}\nu }}))^{{n_{1}}}}\\ {} & ={\sum \limits_{k=0}^{{n_{1}}}}\left(\genfrac{}{}{0.0pt}{}{{n_{1}}}{k}\right){(ab+c)^{k}}{e^{2\pi \mathrm{i}k\nu }}{a^{{n_{1}}-k}}{(1-b{e^{2\pi \mathrm{i}\nu }})^{{n_{1}}-k}}.\end{aligned}\]
Finally, with use of (46) and (45),
\[\begin{aligned}{}& {\int _{-1/2}^{1/2}}\frac{{e^{-2\pi \mathrm{i}{n_{2}}\nu }}{(a+c{e^{2\pi \mathrm{i}\nu }})^{{n_{1}}}}}{{(1-b{e^{2\pi \mathrm{i}\nu }})^{{n_{1}}+1}}}\hspace{0.1667em}d\nu \\ {} & ={\underset{-1/2}{\overset{1/2}{\int }}}\hspace{-0.1667em}\hspace{-0.1667em}\frac{{e^{-2\pi \mathrm{i}{n_{2}}\nu }}}{{(1-b{e^{2\pi \mathrm{i}\nu }})^{{n_{1}}+1}}}{\sum \limits_{k=0}^{{n_{1}}}}\hspace{-0.1667em}\left(\genfrac{}{}{0.0pt}{}{{n_{1}}}{k}\right){(ab+c)^{k}}{e^{2\pi \mathrm{i}k\nu }}{a^{{n_{1}}-k}}{(1-b{e^{2\pi \mathrm{i}\nu }})^{{n_{1}}-k}}d\nu \\ {} & ={\sum \limits_{k=0}^{{n_{1}}}}\left(\genfrac{}{}{0.0pt}{}{{n_{1}}}{k}\right){a^{{n_{1}}-k}}{(ab+c)^{k}}{\int _{-1/2}^{1/2}}\frac{{e^{2\pi \mathrm{i}(k-{n_{2}})\nu }}}{{(1-b{e^{2\pi \mathrm{i}\nu }})^{1+k}}}\hspace{0.1667em}d\nu \\ {} & ={\sum \limits_{k=0}^{\min ({n_{1}},\hspace{0.1667em}{n_{2}})}}\left(\genfrac{}{}{0.0pt}{}{{n_{1}}}{k}\right){a^{{n_{1}}-k}}{(ab+c)^{k}}\left(\genfrac{}{}{0.0pt}{}{{n_{2}}}{k}\right){b^{{n_{2}}-k}}.\end{aligned}\]
 □

Acknowledgement

The author thanks to the referees and editors for careful reading of the paper and providing relevant additional references.

References

[1] 
Adu, N., Richardson, G.: Unit roots test: Spatial model with long memory errors. Statist. Probab. Lett. 140, 126–131 (2018). MR3812259. https://doi.org/10.1016/j.spl.2018.05.003
[2] 
Baran, S.: On the variances of a spatial unit root model. Lith. Math. J. 51(2), 122–140 (2011). MR2805732. https://doi.org/10.1007/s10986-011-9113-9
[3] 
Baran, S., Pap, G.: Parameter estimation in a spatial unilateral unit root autoregressive model. J. Multivariate Anal. 107, 282–305 (2012). MR2890448. https://doi.org/10.1016/j.jmva.2012.01.022
[4] 
Baran, S., Pap, G., Van Zuijlen, M.C.A.: Asymptotic inference for an unstable spatial AR model. Statistics 38(6), 465–482 (2004). MR2109629. https://doi.org/10.1080/02331880412331319297
[5] 
Basu, S., Reinsel, G.C.: Properties of the spatial unilateral first-order ARMA model. Adv. in Appl. Probab. 25(3), 631–648 (1993). MR1234300. https://doi.org/10.2307/1427527
[6] 
Champagnat, F., Idier, J.: On the correlation structure of unilateral AR processes on the plane. Adv. in Appl. Probab. 32(2), 408–425 (2000). MR1778572. https://doi.org/10.1239/aap/1013540171
[7] 
Cressie, N.A.C.: Statistics for Spatial Data. Wiley, New York (1993). MR1239641. https://doi.org/10.1002/9781119115151
[8] 
Jain, A.K.: Advances in mathematical models for image processing. Proc. IEEE 69(5), 502–528 (1981). https://doi.org/10.1109/PROC.1981.12021
[9] 
Kallianpur, G.: Some remarks on the purely nondeterministic property of second order random fields. In: Arató, M., Vermes, D., Balakrishnan, A.V. (eds.) Stochastic Differential Systems, pp. 98–109. Springer, Berlin (1981). https://doi.org/10.1007/BFb0006413 MR0653652
[10] 
Martin, R.J.: A subclass of lattice processes applied to a problem in planar sampling. Biometrika 66(2), 209–217 (1979). MR0548186. https://doi.org/10.1093/biomet/66.2.209
[11] 
Martin, R.J.: The use of time-series models and methods in the analysis of agricultural field trials. Comm. Statist. Theory Methods 19(1), 55–81 (1990). MR1060398. https://doi.org/10.1080/03610929008830187
[12] 
Martin, R.J.: Some results on unilateral ARMA lattice processes. J. Statist. Plann. Inference 50(3), 395–411 (1996). MR1394140. https://doi.org/10.1016/0378-3758(95)00066-6
[13] 
Paulauskas, V.: On unit root for spatial autoregressive models. J. Multivariate Anal. 98(1), 209–226 (2007). MR2292924. https://doi.org/10.1016/j.jmva.2006.08.001
[14] 
Pickard, D.K.: Unilateral Markov fields. Adv. in Appl. Probab. 12(3), 655–671 (1980). MR0578842. https://doi.org/10.2307/1426425
[15] 
Shumway, R.H., Stoffer, D.S.: Time Series Analysis and Its Applications, 4th edn. Springer Texts in Statistics. Springer, Cham (2017). MR3642322. https://doi.org/10.1007/978-3-319-52452-8
[16] 
Tjøstheim, D.: Statistical spatial series modelling. Adv. in Appl. Probab. 10(1), 130–154 (1978). MR0471224. https://doi.org/10.2307/1426722
[17] 
Tjøstheim, D.: Statistical spatial series modelling II: Some further results on unilateral lattice processes. Adv. in Appl. Probab. 15(3), 562–584 (1983). MR0706617. https://doi.org/10.2307/1426619
[18] 
Tjøstheim, D.: Correction: Statistical spatial series modelling. Adv. in Appl. Probab. 16(1), 220 (1984). MR0732140. https://doi.org/10.2307/1427234
[19] 
Tory, E.M., Pickard, D.K.: Unilateral Gaussian fields. Adv. in Appl. Probab. 24(1), 95–112 (1992). MR1146521. https://doi.org/10.2307/1427731
[20] 
Unwin, D.J., Hepple, L.W.: The statistical analysis of spatial series. J. R. Stat. Soc. D 23(3/4), 211–227 (1974). https://doi.org/10.2307/2987581
[21] 
Whittle, P.: On stationary processes in the plane. Biometrika 41(3-4), 434–449 (1954). MR0067450. https://doi.org/10.1093/biomet/41.3-4.434
Reading mode PDF XML

Table of contents
  • 1 Introduction
  • 2 Existence of a stationary solution
  • 3 The autocovariance function
  • 4 Stability and causality
  • 5 Symmetry, nonidentifiability, and special cases
  • 6 Pure nondeterminism
  • 7 Conclusion
  • A Appendix: Auxiliary results
  • Acknowledgement
  • References

Copyright
© 2025 The Author(s). Published by VTeX
by logo by logo
Open access article under the CC BY license.

Keywords
autoregressive models causality discrete random fields purely nondeterministic random fields stationary random fields 60G60 62M10

Metrics
since March 2018
314

Article info
views

165

Full article
views

174

PDF
downloads

41

XML
downloads

Export citation

Copy and paste formatted citation
Placeholder

Download citation in file


Share


RSS

  • Tables
    4
Table 1.
The autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$ for $a=-0.1$, $b=0.5$, $c=0.2$, and ${\sigma _{\epsilon }^{2}}=0.72$
Table 2.
Correspondence of terms of the recurrence equation for fields ${X_{i,j}}$, ${X_{-i,-j}}$, ${X_{-i,j}}$ and ${X_{i,-j}}$
Table 3.
Four sets of parameters $\theta =(a,b,c,{\sigma _{\epsilon }^{2}})$ that determine the same autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$
Table 4.
Relation between the signs of ${f_{1}},\dots ,{f_{4}}$ and line in Table 2 where parameters satisfy the causality condition
Table 1.
The autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$ for $a=-0.1$, $b=0.5$, $c=0.2$, and ${\sigma _{\epsilon }^{2}}=0.72$
Values of ${\boldsymbol{\gamma }_{\boldsymbol{X}}}\boldsymbol{(}{\boldsymbol{h}_{\mathbf{1}}}\mathbf{,}{\boldsymbol{h}_{\mathbf{2}}}\boldsymbol{)}$
${\boldsymbol{h}_{\mathbf{2}}}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\boldsymbol{-}\mathbf{2}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\boldsymbol{-}\mathbf{1}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\mathbf{0}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\mathbf{1}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\mathbf{2}$ ${\boldsymbol{h}_{\mathbf{1}}}\mathbf{=}\mathbf{3}$
3 0 0 0.125 0.1125 0.0225 –0.00225
2 0 0 0.25 0.15 0.0075 –0.003
1 0 0 0.5 0.15 –0.015 0.0015
0 0 0 1 0 0 0
–1 –0.015 0.15 0.5 0 0 0
–2 0.0075 0.15 0.25 0 0 0
–3 0.0225 0.1125 0.125 0 0 0
Table 2.
Correspondence of terms of the recurrence equation for fields ${X_{i,j}}$, ${X_{-i,-j}}$, ${X_{-i,j}}$ and ${X_{i,-j}}$
Line Field of Coefficients Error field
no., $\boldsymbol{m}$ interest variance values
1 ${X_{i,j}}$ a b c ${\sigma _{\epsilon }^{2}}$ ${\epsilon _{i,j}}$
2 ${X_{-i,-j}}$ $-b/c$ $-a/c$ $1/c$ ${c^{-2}}{\sigma _{\epsilon }^{2}}$ $-{c^{-1}}{\epsilon _{1-i,1-j}}$
3 ${X_{-i,j}}$ $1/a$ $-c/a$ $-b/a$ ${a^{-2}}{\sigma _{\epsilon }^{2}}$ $-{a^{-1}}{\epsilon _{1-i,j}}$
4 ${X_{i,j-1}}$ $-c/b$ $1/b$ $-a/b$ ${b^{-2}}{\sigma _{\epsilon }^{2}}$ $-{b^{-1}}{\epsilon _{i,1-j}}$
Table 3.
Four sets of parameters $\theta =(a,b,c,{\sigma _{\epsilon }^{2}})$ that determine the same autocovariance function ${\gamma _{X}}({h_{1}},{h_{2}})$
a b $c=-ab$ ${\sigma _{\epsilon }^{2}}$
${a^{-1}}$ ${b^{-1}}$ $-{a^{-1}}{b^{-1}}$ ${a^{-2}}{b^{-2}}{\sigma _{\epsilon }^{2}}$
${a^{-1}}$ b $-{a^{-1}}b$ ${a^{-2}}{\sigma _{\epsilon }^{2}}$
a ${b^{-1}}$ $-a{b^{-1}}$ ${b^{-2}}{\sigma _{\epsilon }^{2}}$
Table 4.
Relation between the signs of ${f_{1}},\dots ,{f_{4}}$ and line in Table 2 where parameters satisfy the causality condition
Signs of ${\boldsymbol{f}_{\mathbf{1}}}\mathbf{,}\dots \mathbf{,}{\boldsymbol{f}_{\mathbf{4}}}$ The Which Outcome
signs parameter-
imply ization
that is causal, $\boldsymbol{m}$
${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$, ${f_{4}}\gt 0$ 1 ${X_{i,j}}$ is causal w.r.t. ${\epsilon _{i,j}}$
${f_{1}}\lt 0$, ${f_{2}}\gt 0$, ${f_{3}}\gt 0$, ${f_{4}}\lt 0$ $c\gt 0$ 2 ${X_{-i,-j}}$ is causal
${f_{1}}\gt 0$, ${f_{2}}\lt 0$, ${f_{3}}\lt 0$, ${f_{4}}\gt 0$ $c\lt 0$ w.r.t. ${\epsilon _{1-i,1-j}}$
${f_{1}}\lt 0$, ${f_{2}}\lt 0$, ${f_{3}}\gt 0$, ${f_{4}}\gt 0$ $a\gt 0$ 3 ${X_{-i,j}}$ is causal
${f_{1}}\gt 0$, ${f_{2}}\gt 0$, ${f_{3}}\lt 0$, ${f_{4}}\lt 0$ $a\lt 0$ w.r.t. ${\epsilon _{1-i,j}}$
${f_{1}}\lt 0$, ${f_{2}}\gt 0$, ${f_{3}}\lt 0$, ${f_{4}}\gt 0$ $b\gt 0$ 4 ${X_{i,-j}}$ is causal
${f_{1}}\gt 0$, ${f_{2}}\lt 0$, ${f_{3}}\gt 0$, ${f_{4}}\lt 0$ $b\lt 0$ w.r.t. ${\epsilon _{i,1-j}}$

MSTA

MSTA

  • Online ISSN: 2351-6054
  • Print ISSN: 2351-6046
  • Copyright © 2018 VTeX

About

  • About journal
  • Indexed in
  • Editors-in-Chief

For contributors

  • Submit
  • OA Policy
  • Become a Peer-reviewer

Contact us

  • ejournals-vmsta@vtex.lt
  • Mokslininkų 2A
  • LT-08412 Vilnius
  • Lithuania
Powered by PubliMill  •  Privacy policy