# Yule-Walker Estimation for the Moving-Average Model

Authors: Chrysoula Dimitriou-Fakalou
Journal: International Journal of Stochastic Analysis
 (2011)

---

## Abstract

The standard Yule-Walker equations, as they are known for an autoregression, are generalized
to involve the moments of a moving-average process indexed on any number of dimensions. Once
observations become available, new moments estimators are set to imitate the theoretical equations.
These estimators are not only consistent but also asymptotically normal for any number of indexes.
Their variance matrix resembles a standard result from maximum Gaussian likelihood estimation.
A simulation study is added to conclude on their efficiency.


---

## Body

## 1. Introduction

At first, processes taking place on more than one dimension were considered on surfaces or in the three-dimensional space for the sake of spatial dynamics. Next, the time series Autoregressive Moving-Average model could be as well used to clothe the covariance dependence of spatial point processes. Nevertheless, the assumptions of causality and invertibility were based on conventional orderings of the vector indexes, as these were introduced by Whittle [1] and generalized by Guyon [2], and they failed to look natural on the spatial axes. The more general bilateral equations considered by Whittle [1] also related to serious problems during the estimation of the parameters. Besag’s [3] automodels as well as Besag’s [4] equations for observations on irregular sites moved in a new direction but the weaknesses during the estimation remained.The inclusion of the time axis has provided more convenience for the spatial analysis that may be summarized in two ways. On the one hand when the observations are collected on irregular sites, the asymptotic properties of the estimators for the spatial dependence parameters can now be established as the number of timings only increases to infinity; for example, Dimitriou-Fakalou [5] has proposed a methodology on how to proceed. On the other hand when the regular recordings are increasing towards all the spatial as well as the time axes, the unilateral ordering of indexes can now be arranged in a meaningful way and the causal and invertible ARMA model may be preferred as a natural spatiotemporal equation.Special cases of the ARMA model are the pure autoregressive and moving-average equations and they both relate to useful properties as in the spectral density of interest. The autoregressive process uses a spectral density with “finite denominator”, which for Gaussian processes [3] translates into a finite conditional expectation of any value based on all other values. A moving-average process uses a “finite numerator” in the spectral density, resulting in the autocovariance function with nonzero values on a finite set of “lags” only.For a causal and finite autoregressive equation, the method of moments estimators according to the Yule-Walker equations almost coincides with the least squares or maximum Gaussian likelihood estimators and their consistency as well as their asymptotic normality may be established. Nevertheless, if the estimation for a pure moving-average or ARMA multi-indexed process takes place using its infinite autoregressive representation, a problem known as the “edge-effect” [2] resurrects and the asymptotic normality cannot be guaranteed. A modification of the Gaussian likelihood, which has been written from the autoregressive representation of each observation based on its “previous” ones, has been proposed by Yao and Brockwell [6] to confine the edge-effect with a weak condition of a finite second moment, but this works for the case of processes indexed on two dimensions only.In this paper, for the multidimensional processes that may be clothed via an invertible and finite moving-average equation, a new method of estimation is proposed imitating the Yule-Walker equations. The finite number of nonzero autocovariances of the moving-average process is used most advantageously to eliminate the edge-effect. As a result, the proposed estimators are not only consistent but also asymptotically normal using relaxed conditions of finite second moments.


## 2. Theoretical Yule-Walker Equations

Forℤ={0,±1,±2,…} and d any positive integer, real-valued processes indexed on ℤd will be considered, such as time series or spatial processes on one line transect (d=1), spatial processes in the plane (d=2) or in the three-dimensional space, as well as the combination of spatiotemporal processes on ℤ2, ℤ3, or ℤ4, respectively. For any v∈ℤd, an invertible moving-average process {Y(v)} will be defined by the following equation:(2.1)Y(v)=ε(v)+∑n=1qθinε(v-in),
where {ε(v)} are uncorrelated, zero-mean random variables with positive and finite variance σ2. In (2.1) a unilateral ordering of the fixed lags 0<i1<⋯<iq has been considered. For d=1, this is the standard ordering of integer numbers and, for d=2, the conventional unilateral ordering has been introduced by Whittle [1]. Guyon [2] has generalized the ordering of two vector indexes on ℤd, d≥3. As for the invertibility condition, it is only imposed to make sure that it can be written ε(v)=Y(v)+∑i>0ΘiY(v-i), ∑i>0|Θi|<∞ also using the ordering i>0. Anderson and Jury  [7] have provided a sufficient condition on the invertibility of such filters.For the complex numbersz1,…,zd and the vector z=(z1,…,zd), the polynomial  θ(z)≡1+∑n=1qθinzin is considered, where zi=z(i1,…,id)≡z1i1,…,zdid for i1,…,id∈ℤ. Similarly, the backwards operators B1,…,Bd and the vector backwards operator B=(B1,…,Bd) are considered, such that Biε(v)≡ε(v-i). So (2.1) may now be rewritten as(2.2)Y(v)=θ(B)ε(v).
Writing z-1≡(z1-1,…,zd-1), the autoregressive sequence {X(v)} may be defined from {ε(v)} according to the following:(2.3)θ(B-1)X(v)≡X(v)+∑n=1qθinX(v+in)=ε(v).
The auto-regression is unilateral but not causal, in the sense that it can be written that X(v)=ε(v)+∑i>0Θiε(v+i), ∑i>0|Θi|<∞ as a function of the “future” terms ε(v+i), i≥0. Moreover, writing the two polynomials(2.4)γ(z)=θ(z)θ(z-1)≡∑i∈Fγizi,c(z)=γ(z)-1≡∑i∈Zdcizi,
it holds that(2.5)Y(v)=γ(B)X(v)orX(v)=c(B)Y(v).
The set ℱ≡{in,-in,in-im,n,m=1,…,q} in (2.4) has members the vector “lags”, where {Y(v)} have nonzero autocorrelations. The set {in,in-im,n,m=1,…,q,n>m} that uses the “positive” lags of ℱ will have cardinality q*.The relationship in (2.5) provides the way to write Yule-Walker equations for a moving-average process. The processes {Y(v)} and {X(v)} are not only “caused” by the same sequence {ε(v)}, one from its “past” and the other from its “future” values, but also have “reverse” second-order dependence; that is, the spectral density of the auto-regression is(2.6)gX(ω1,…,ωd)≡σ2(2π)d1θ(z)θ(z-1)=σ2(2π)dc(z),zk=e-iωk,ωk∈(-π,π),k=1,…,d,
where according to (2.4), E(X(v)X(v-i))=σ2ci and the spectral density of the moving-average process is(2.7)gY(ω1,…,ωd)≡σ2(2π)dθ(z)θ(z-1)=σ2(2π)dγ(z),zk=e-iωk,ωk∈(-π,π),k=1,…,d,
which also implies that E(Y(v)Y(v-i))=σ2γi. While the original Yule-Walker equations for {X(v)} dictate(2.8)c0+∑n=1qθincin=1,(2.9)ci+∑n=1qθinci-in=0,i>0,
from (2.9), the following (q+1) equations may be written for i>0, that is, (2.10)ci+θi1ci-i1+θi2ci-i2+⋯+θiqci-iq=0,θi1ci+i1+θi12ci+θi1θi2ci+i1-i2+⋯+θi1θiqci+i1-iq=0,⋮θiqci+iq+θiqθi1ci+iq-i1+θiqθi2ci+iq-i2+⋯+θiq2ci=0
and their sum ∑j∈ℱγjci-j=∑j∈ℱγjcj-i=∑-j∈ℱγ-jc-j-i=∑j∈ℱγjcj+i=0, which may be rewritten as(2.11)∑j+i∈Fcjγj+i=0,i≠0.
Further to (2.11), it holds that(2.12)∑j∈Fcjγj=1
from (2.8) and (2.9). It is true that (2.11) and (2.12) are equivalent to γ(z)c(z)=1. Like the standard Yule-Walker equations, the general versions (2.11) and (2.12) are constructed from the coefficients θi1,…,θiq. As the finite number of nonzero autocovariances of {Y(v)} is involved now, they will be used as theoretical prototypes to be imitated by their sample analogues for the sake of estimation of the moving-average parameters.Finally, as the assumption of invertibility of the moving-average equation used for parameterization will be essential, it should be noted that Dimitriou-Fakalou [8] has presented in detail the merits of a causal and invertible spatiotemporal ARMA model. As a special case, the moving-average process {Yt(s)}, where t is the time index and s are spatial indexes, may satisfy the equation Yt(s)=εt(s)+∑i=1qt∑j∈ℳiθi,jεt-i(s-j), where {εt(s)} are zero mean, uncorrelated variables with identical, finite variance, qt is a finite positive integer, ℳ1,…,ℳqt have a finite cardinality, and it holds that(2.13)εt(s)=Yt(s)+∑i>0∑jΘi,jYt-i(s-j),∑i>0∑j|Θi,j|<∞;
that is, the model is invertible. Note that, in either equation, no unnatural ordering is forced over space as the unilateral conventions “run” over the past timings only. Hence the model both uses a space-time unilateral ordering of indexes and is meaningful. As a result, for the special case of a stationary process with autocovariance function identified to take zero values after a reasonably small number of lags on all spatial and temporal dimensions, the invertible moving-average equation should be preferred and the new results presented next will be useful.


## 3. Method of Moments Estimators

For𝒮⊂ℤd being a set of finite cardinality N, there are available {Y(v),v∈𝒮} from (2.2) with true parameters θ0=(θi1,0,…,θiq,0)τ that are to be estimated. The following condition will be assumed to be true.Condition C.
The parameter spaceΘ⊂ℝq is a compact set containing the true value θ0 as an inner point. Further, for any θ=(θi1,…,θiq)τ∈Θ the moving-average model (2.2) is invertible.Based on𝒮, a set ℱv will be defined for any v∈ℤd. More specifically, it will hold that i∈ℱv only if v-i∈𝒮. Next, for the corrected set 𝒮* with N* elements, it will hold v∈𝒮* if v+i∈𝒮 for every i∈ℱ. Thus for any v∈𝒮*, it holds that ℱ⊆ℱv. The reduction from 𝒮 to 𝒮* is essential as it will guarantee later that the edge-effect will fade away; this is because the set 𝒮*⊂𝒮 includes the locations v that have all their “neighbors” v-i with E(Y(v)Y(v-i))≠0 being available. As the source of the edge-effect is the speed of the absolute bias to zero, 𝒮* will guarantee that whatever “is missing” from the finite sample will not add at all to the bias of the estimators.By imitating (2.11) the estimators θ̂=(θ̂i1,…,θ̂iq)τ are set to be the solutions of equations:(3.1)∑v∈S*{∑j-in∈FvĉjY(v+in-j)}Y(v)≡0,n=1,…,q,
where θ̂(z)≡1+∑n=1qθ̂inzin and ĉ(z)≡∑i∈ℤdĉizi≡{θ̂(z)θ̂(z-1)}-1. Using Proposition 3.1, the equations(3.2)∑v∈S*{∑j-(in-im)∈FvĉjY(v+in-im-j)}Y(v)≡0,n,m=1,…,q,n>m,
will be used additionally, in order to make sure that the estimators are consistent. The proposition guarantees that the use of q⋆ equations used as prototypes for estimation can provide a unique way to compute the q parameters of the invertible equation.Proposition 3.1.
The polynomialsθ0(z)=1+∑n=1qθin,0zin and θ0(z)-1=1+∑i>0Θi,0zi with ∑i>0|Θi,0|<∞ as well as γ0(z)=θ0(z)θ0(z-1)≡∑j∈ℱγj,0zj and c0(z)=γ0(z)-1≡∑j∈ℤdcj,0zj are considered. For any θ=(θi1,…,θiq)τ∈Θ, such that it holds that θ(z)=1+∑n=1qθinzin and θ(z)-1=1+∑i>0Θizi with ∑i>0|Θi|<∞, as well as γ(z)=θ(z)θ(z-1)≡∑j∈ℱγjzj and c(z)=γ(z)-1≡∑j∈ℤdcjzj, the unique solution that satisfies the q⋆ equations ∑j∈ℱγj,0cj-in=0, n=1,…,q, and ∑j∈ℱγj,0cj-(in-im)=0, n,m=1,…,q, n>m, is θin≡θin,0, n=1,…,q.For mathematical convenience, new variables depending on the sampling set are definedHY(v)≡Y(v), v∈𝒮 and HY(v)≡0, v∉𝒮, and (3.1) may be rewritten as(3.3)∑v∈S*{∑j∈Zdcj,0HY(v+in-j)}Y(v)-Jn(θ̂-θ0)=0,n=1,…,q,
or(3.4)∑v∈S*{c0(B)HY(v+in)}Y(v)-Jn(θ̂-θ0)=0,n=1,…,q,
where the zero subindex relates to θ0. In (3.3) and (3.4), it is Jn=(Jn,1,…,Jn,q) with elements(3.5)Jn,m=∑v∈S*{c0(B)(θ0(B)-1HY(v+in-im)+θ0(B-1)-1HY(v+in+im))}Y(v)+OP(N∥θ̂-θ0∥2),n,m=1,…,q.
The last term in (3.5) will become more obvious when the properties of the estimators will be established but for now it might be justified by a second-order Taylor expansion, where it holds(3.6)(θ̂in-θin,0)(θ̂im-θim,0)‖θ̂-θ0‖2⟶λn,m∈(0,∞),n,m=1,…,q
due to the fact that θ̂ will be shown to be consistent estimators in Theorem 4.1. Further, the sums over 𝒮* involving the second derivatives of the Taylor expansion, when divided by N* or N (N*/N→1), will also tend in probability to the relevant finite expectations when {ε(v)} will be assumed to be independent and identically distributed random variables with a finite variance.Finally by imitating (2.12), the estimator of the error variance may be defined as(3.7)σ̂2≡∑j∈Fĉj∑v∈S*Y(v)Y(v-j)N*.Example 3.2.
For the simplest one-dimensional case whenq*=q=1, it is demonstrated in detail how to compute the new estimator; some other cases will briefly be considered later. Suppose that are available {Y(1),…,Y(N)} from Y(v)=e(v)+θe(v-1), where |θ|<1, θ≠0, and {e(v)} are uncorrelated random variables with variance unity. The original set 𝒮≡{1,2,…,N} will be corrected to 𝒮*≡{2,3,…,N-1} that it is the maximal set, such that for all v∈𝒮* it holds that E(Y(v)Y(v*))=0, v*∉𝒮. Note that the “lags” on which the autocorrelations of Y are not equal to zero are as in the set ℱ={1,0,-1}; further, for any v∈ℤ it holds that ℱv={v-1,…,v-N} and, thus, for v∈𝒮* it is true that ℱ⊆ℱv.
The autocovariance of the relevant auto-regression at the “lag”j∈ℤ is (-θ)|j|/(1-θ2). As a result, the estimator θ̂ will be the solution of the following equation:
(3.8)0=-θ̂∑v=2N-1Y(v)2+∑v=2N-1Y(v)Y(v+1)+θ̂2∑v=2N-1Y(v)Y(v-1)-θ̂∑v=2N-2Y(v)Y(v+2)-θ̂3∑v=3N-1Y(v)Y(v-2)+θ̂2∑v=2N-3Y(v)Y(v+3)+θ̂4∑v=4N-2Y(v)Y(v-3)-⋯+(-θ̂)N-4Y(2)Y(N-1)+(-θ̂)N-2Y(N-2)Y(1).
Furthermore, when {e(v)} are independent and identically distributed random variables, it will hold in probability that ∑vY(v)Y(v-i)/N→0 for |i|=2,3,… and the estimator may be approximated as a solution of the quadratic equation:
(3.9)0≡-θ̂∑v=2N-1Y(v)2+∑v=2N-1Y(v)Y(v+1)+θ̂2∑v=2N-1Y(v)Y(v-1).
Later, it will be investigated whether the approximation is worthwhile. The estimator reduces to
(3.10)θ̂=∑v=2N-1Y(v)2±(∑v=2N-1Y(v)2)2-4(∑v=2N-1Y(v)Y(v-1))(∑v=2N-1Y(v)Y(v+1))2∑v=2N-1Y(v)Y(v-1)
or
(3.11)θ̂=1±1-4ρ̂+ρ̂-2ρ̂+
with
(3.12)ρ̂+=∑v=2N-1Y(v)Y(v-1)∑v=2N-1Y(v)2,ρ̂-=∑v=2N-1Y(v)Y(v+1)∑v=2N-1Y(v)2.
For the actual nonzero auto-correlation ρ=θ/(1+θ2), it holds that |ρ|<1/2 and D=1-4ρρ>0. As a result, if 0<θ<1 and 0<ρ<1/2, then 1/2ρ>1, and the value
(3.13)1+1-4ρρ2ρ=12ρ+D2ρ
is bigger than 1. If, on the other hand, -1<θ<0, and -1/2<ρ<0 then 1/2ρ<-1 and the same value is smaller than −1. Thus, the estimator is approximately
(3.14)θ̂=∑v=2N-1Y(v)2-(∑v=2N-1Y(v)2)2-4(∑v=2N-1Y(v)Y(v-1))(∑v=2N-1Y(v)Y(v+1))2∑v=2N-1Y(v)Y(v-1).
Although the distribution of θ̂ is not known for small sample sizes and |ρ̂+|,|ρ̂-| risk to be outside the (-0.5,0.5) range with some positive probability, both the consistency and the asymptotic normality will be established next. As opposed to what is known for the maximum Gaussian likelihood estimators, the proposed estimators will be asymptotically normal even for processes on more than one dimension. Thus, the new methods might be more useful for higher number of dimensions as well as large sample sizes.


## 4. Properties of Estimators

The asymptotic normality of standard estimators for the parameters of stationary processes indexed onℤd, d≥2, has not been established yet. More specifically as Guyon [2] demonstrated, the bias of the genuine Gaussian likelihood estimators, computed from N regular recordings increasing to infinity on all dimensions and at equal speeds, is of order N-1/d. Nevertheless in order to secure the asymptotic normality, the absolute bias when multiplied by N1/2 should tend to zero, which is only true for d=1 as in Brockwell and Davis [9].Regarding the proposed estimators of this paper, their bias will stem from (3.3) as in the expected value of ∑v∈𝒮*{∑j∈ℤdcj,0(Y(v+in-j)-HY(v+in-j))Y(v)}/N*, n=1,…,q, which expresses what is “missing” from the sample. Even when the bias is multiplied by N*1/2, the random variable N*-1/2∑v∈𝒮*∑j-in∉ℱvcj,0Y(v+in-j)Y(v), n=1,…,q has a zero mean due to the selection 𝒮* and the fact that the “finite” autocovariance function of a moving-average process is being used. Thus, the edge-effect will not be an obstacle to establish the asymptotic normality of the estimators.A series of ARMA model-based arguments have been used before to deal with the edge-effect by Yao and Brockwell [6] for two-dimensional processes as well as a weak condition of finite second moments. Guyon [2], on the other hand, used the periodogram as in the Gaussian likelihood of observations proposed by Whittle [1] and required that the fourth moments would be finite. Dahlhaus and Künsch [10] improved certain weaknesses of Guyon's [2] modification on the likelihood but paid the price of having to reduce the dimensionality to secure their results. The proposed estimators of this paper will also use the weak condition of finite second moments, which is a great advantage, since they will refer to moving-average processes on any number of dimensions.Condition 4.
(i) For a set𝒮≡𝒮N⊂ℤd of cardinality N, it is written N→∞ if the length M of the minimal hypercube including 𝒮, say 𝒮⊆𝒞M, and the length m of the maximal hypercube included in 𝒮, say 𝒞m⊆𝒮, are such that M,m→∞. (ii) Further, as M,m→∞ it holds that M/m is bounded away from ∞.Theorem 4.1.
If{ε(v)}~IID(0,σ2), then under Conditions C1 and C2(i) and as N→∞, it holds that
(4.1)θ̂⟶Pθ0,σ̂2⟶Pσ2.
From (3.3) and (3.4), all the q equations together may be written as
(4.2)N1/2(θ̂-θ0)=(JN)-1(N-1/2∑v∈S*HY(v)),
where Jτ=(J1τ,…,Jqτ) and
(4.3)HY(v)≡(c0(B)HY(v+i1)⋮c0(B)HY(v+iq))Y(v),v∈Zd,
which depends on the sampling set 𝒮. The next proposition reveals what happens to the part (J/N) in (4.2). Then Theorem 4.3 establishes the asymptotic normality of the estimators.Proposition 4.2.
Let the polynomialθ0(z)-1=(1+∑n=1qθin,0zin)-1≡∑i≥0Θi,0zi,Θ0,0. If {ε(v)}~IID(0,σ2), then under Conditions C1 and C2(i) and as N→∞, it holds that
(4.4)JN→Pσ2Θ0≡σ2(100⋯0Θi2-i1,010⋯0⋮⋱Θiq-i1,0Θiq-i2,0Θiq-i3,01).Theorem 4.3.
For{W(v)}~IID(0,1), the auto-regression {η(v)} is defined by θ0(B)η(v)≡W(v). Also the vector ξ≡(η(-i1),…,η(-iq))τ and the variance matrix
(4.5)Wq*≡Var⁡(ξ∣W(-i1-i,i>0,i≠i2-i1,…,iq-i1))
are considered. If {ε(v)}~IID(0,σ2), then under Conditions C1 and C2 and as N→∞, it holds that
(4.6)N1/2(θ̂-θ0)→DN(0,Wq*-1).


## 5. Empirical Results

In this section, an empirical comparison of the proposed estimators to maximum Gaussian likelihood estimators for moving-average processes is presented in detail. The theoretical foundations in favour of the new estimators have been provided already when their asymptotic normality on any number of dimensions has been established based on finite second moments only. As a result, the speed of the bias to zero is not expected to cause them to perform worse than maximum likelihood estimators, especially when the dimensionality is large.On the other hand, Theorem4.3 attributes to the new estimators the variance matrix Wq*-1 when, according to Hannan [11], efficient estimation results in Wq-1 with Wq≡Var(ξ) and the same notation as in Theorem 4.3. Wq and Wq* are defined similarly but they are not, in general, equal. It seems that as the number of moving-average parameters of the process increases, the two types elements get closer. Nevertheless, as the pure moving-average model will be preferred when it is parsimonious, a decision is made here whether the new estimators are efficient then.The investigation has started with the one-dimensional case by generating{Y(1),…,Y(N)} from the model with one parameter only. The moments estimator has been approximated as in the example earlier, while the minimizer of ∑v∈𝒮*e*2(v) with respect to the parameter θ has been considered to be the Gaussian likelihood estimator. The true values θ=0.5,-0.8 have been considered and the sample size has been set equal to N=30,100,200,300. Very encouraging results for the efficiency of the proposed estimator are presented in Table 1 as even when the sample size is still small, extreme differences in the variances of the two types estimators cannot be detected. It is the bias of the moments estimator that seems to be the only reason why it might be outperformed by the likelihood estimator in terms of the Mean Squared Error. Nevertheless, the speed with which the bias tends to zero is much faster as one would expect and, eventually, the new estimator for θ=0.5 performs better altogether.Table 1
Approximate bias, variance, and Mean Squared Error of the Method of Moments and the Maximum Likelihood Estimator for the MA(1) parameter from 1000 replications.
θ=0.5θ=-0.8N=30BiasVarMSEBiasVarMSEMME−0.4033450.000209720.1628970.09310520.000979120.00964771MLE0.017950.000159840.000482040.0209490.000623470.00104331N=100BiasVarMSEBiasVarMSEMME−0.1238120.000214440.01554390.02179280.000979120.00145404MLE0.005210.000313290.000340430.005990.000639360.00067524N=200BiasVarMSEBiasVarMSEMME−0.04226530.000211730.001998080.01232750.000979120.00113109MLE−0.001160.000220680.000222020.002750.000623480.00063104N=300BiasVarMSEBiasVarMSEMME−0.007088140.000237680.00028792−0.001665970.000979120.00098189MLE0.003130.000291310.000361150.003860.000655440.00067034Next for the two-dimensional case,{Y(u,v),u,v=1,…,n} have been generated from(5.1)Y(u,v)=e(u,v)+βe(u-1,v)+γe(u,v-1),|β|+|γ|<1,
with {e(u,v),u,v∈ℤ} being independent, Gaussian random variables with zero mean and variance unity. Using similar arguments as in the example, the moments estimators of β and γ might be derived by solving simultaneously (5.2)0≡∑Y(u,v)Y(u+1,v)ĉ(0,0)+∑Y2(u,v)ĉ(-1,0)+∑Y(u,v)Y(u-1,v)ĉ(-2,0)+∑Y(u,v)Y(u,v+1)ĉ(-1,1)+∑Y(u,v)Y(u,v-1)ĉ(-1,-1)+∑Y(u,v)Y(u+1,v-1)ĉ(0,-1)+∑Y(u,v)Y(u-1,v+1)ĉ(-2,1),0≡∑Y(u,v)Y(u,v+1)ĉ(0,0)+∑Y2(u,v)ĉ(0,-1)+∑Y(u,v)Y(u,v-1)ĉ(0,-2)+∑Y(u,v)Y(u+1,v)ĉ(1,-1)+∑Y(u,v)Y(u-1,v)ĉ(-1,-1)+∑Y(u,v)Y(u+1,v-1)ĉ(1,-2)+∑Y(u,v)Y(u-1,v+1)ĉ(-1,0),
as well as the third equation (5.3)0≡∑Y(u,v)Y(u+1,v-1)ĉ(0,0)+∑Y2(u,v)ĉ(-1,1)+∑Y(u,v)Y(u-1,v+1)ĉ(-2,2)+∑Y(u,v)Y(u+1,v)ĉ(0,1)+∑Y(u,v)Y(u-1,v)ĉ(-2,1)+∑Y(u,v)Y(u,v+1)ĉ(-1,2)+∑Y(u,v)Y(u,v-1)ĉ(-1,0),
which is needed to secure the consistency (Proposition 3.1). All of the sums above extend over 𝒮*≡{2,…,n-1}×{2,…,n-1}. The estimators introduced here are β̂,γ̂, and in the three equations, it has been considered that ĉ(j,k)≡Corr(X(u,v),X(u-j,v-k)), where(5.4)X(u,v)+β̂X(u-1,v)+γ̂X(u,v-1)=e*(u,v)
and {e*(u,v)} are zero mean, uncorrelated random variables all with variance unity. Under the condition |β̂|+|γ̂|<1, a special form of “causality” X(u,v)≡∑j,k=0∞Θ̂(j,k)e*(u-j,v-k),∑j,k=0∞|Θ̂(j,k)|<∞ is implied and the standard Yule-Walker equations may be used to easily write the needed autocorrelations ĉ, as functions of β̂ and γ̂. Again, the maximum likelihood estimators have been approximated by the minimizers of ∑u,v(Y(u,v)-βe(u-1,v)-γe(u,v-1))2 with respect to β and γ.Table2 verifies once more the conclusions drawn from the one-dimensional case. It is safe to consider that the new estimators are efficient and this is very apparent in the case that the parameters are in absolute value further away from zero (β=0.5 and γ=0.45). The striking case with the small sample size of around 50 points observed in the plane reveals that the variances of the moments estimators may even happen to be smaller. On the other hand, the bias heavily affects the results for the MSE for small sample sizes. Nevertheless as the sample size increases, the absolute bias of the likelihood estimators does not seem to decrease at all versus the bias of the proposed estimators that speedily reaches the zero value. As a result, the new estimators eventually equalize the MSE performance of the standard estimators.Table 2
Approximate bias, variance, and mean squared error of the method of moments and the maximum likelihood estimators (based onN=n2 observations) for the two-dimensional MA parameters β and γ from 1000 replications.
β=0.5γ=0.45n=7BiasVarMSEBiasVarMSEMME0.028420.000108790.00091649−0.078580.000012090.0061869MLE0.00060.00020230.00020266−0.006040.000062440.00009892n=30BiasVarMSEBiasVarMSEMME0.01490.000371730.00059374−0.007460.00020230.00025795MLE−0.001840.00030220.000305580.000660.000167930.00016837n=100BiasVarMSEBiasVarMSEMME0.005760.000259840.000293020.00650.000239860.00028211MLE−0.000440.000220680.00022087−0.000240.00020230.00020236β=0.15γ=-0.3n=30BiasVarMSEBiasVarMSEMME0.026380.000044060.000739960.032740.000016880.00108879MLE0.000120.000022480.000022490.000040.000108790.00010879n=100BiasVarMSEBiasVarMSEMME−0.04226530.000211730.001998080.01232750.000979120.00113109MLE−0.000580.000022480.000022810.000260.000004020.00008408



---

*Source: 151823-2011-08-14.xml*