One-Step Recurrences for Stationary Random Fields on the Sphere

Recurrences for positive definite functions in terms of the space dimension have been used in several fields of applications. Such recurrences typically relate to properties of the system of special functions characterizing the geometry of the underlying space. In the case of the sphere ${\mathbb S}^{d-1} \subset {\mathbb R}^d$ the (strict) positive definiteness of the zonal function $f(\cos \theta)$ is determined by the signs of the coefficients in the expansion of $f$ in terms of the Gegenbauer polynomials $\{C^\lambda_n\}$, with $\lambda=(d-2)/2$. Recent results show that classical differentiation and integration applied to $f$ have positive definiteness preserving properties in this context. However, in these results the space dimension changes in steps of two. This paper develops operators for zonal functions on the sphere which preserve (strict) positive definiteness while moving up and down in the ladder of dimensions by steps of one. These fractional operators are constructed to act appropriately on the Gegenbauer polynomials $\{C^\lambda_n\}$.


Introduction
This paper develops operators for zonal functions on the sphere which preserve (strict) positive definitenesss while moving up and down in the ladder of dimensions by steps of one. The operators provide tools for forming families of (strictly) positive definite zonal functions. Such (strictly) positive definite zonal functions can be used as covariance models for estimating regionalized variables and also for interpolation on spheres.
Within a deterministic context, zonal positive definite functions on the sphere have been used for interpolation or approximation of scattered data (see [10,11] and the references therein). The standard ansatz in this setting is a linear combination of spherical translates of a fixed (zonal) basis function. While the present paper could well have been stated within the context of approximation on the sphere, we rather chose to provide a probabilistic framework, which is to some extent is equivalent, i.e., the theory of regionalized variables.
Regionalized variables on spherical domains can nicely be modeled using random fields on spheres [6,15]. Such a random field is given through a set of random variables, Z(x) say, where x ∈ S d−1 . Assuming the field to be Gaussian, i.e., for every n ∈ N, (Z(x 1 ), . . . , Z(x n )) T has a multivariate Gaussian distribution for any choice of x 1 , . . . , x n ∈ S d−1 , the distribution can be characterized by its first two moments.
Assuming second order (weak) stationarity, the covariance for an isotropic model is determined by a function Cov Z(x), Z(y) = f (cos θ), x, y ∈ S d−1 , where θ = θ(x, y) = arccos x T y is the geodesic distance between the points x and y on the sphere S d−1 .
As a consequence of Kolmogorov's extension theorem (see [5,Theorem 36.3]), the set of isotropic Gaussian random fields can be identified with the set of zonal positive definite functions on the sphere [12,15]. We note in passing that Lévy named such processes Brownian motion. c i c j g(θ(x i , x j )) ≥ 0 holds true for all c 1 , . . . , c n ∈ R. The function is (zonal) strictly positive definite on S d−1 if the inequality holds in the strict sense for all c 1 , . . . , c n ∈ R not vanishing simultaneously.
Although the natural distance on the unit sphere is an angle in [0, π], it is convenient for the purpose of this paper to consider functions in x = cos θ ∈ [−1, 1], instead. Thus, by Λ d−1 we will denote the cone of all functions f ∈ C[−1, 1] such that f (cos · ) is positive definite on S d−1 . Λ + d−1 will denote the subcone of all strictly positive definite functions in Λ d−1 . Gaussian random fields have been widely applied to statistically analyze spatial phenomena [8,16,23]. In particular, kriging allows prediction of spatial variables from given samples at arbitrary locations. The key ingredient for such an approach lies in determining a suitable model for the covariance function of the spherical random field. Commonly, such a model can be inferred from given data through fitting a parametric family of models (i.e., estimation of the covariance).
Models for covariance functions have further been used for simulation of stationary random fields. Matheron [17] suggested a method based on proper averaging of stationary random fields on a lower dimensional space. In the Euclidean setting this so-called turning bands method works as follows: Given a stationary random field Z 1 on the real line with covariance function C 1 and a randomly chosen direction ξ ∈ S d−1 , Z ξ = Z 1 (x T ξ) defines a stationary random field on R d . Averaging over all directions ξ leads to a stationary field on R d the covariance function of which, C d say, relates to C 1 via the so-called turning bands operator The turning bands operator represents one example out of a suite of operators, mapping radial positive definite functions on R d onto such functions on a higher or lower dimensional space. Wendland [25], Wu [26] and Gneiting [13] used such operators to derive compactly supported functions of a given smoothness. Recurrences for radial positive definite functions in general have been studied by several authors [22,27]. Due to Schoenberg's characterization of radial positive definite functions and the fact that scale mixtures of such functions preserve positive definiteness, recurrence operators can be derived from corresponding relations between special functions. In the case of radial functions on R d , the appropriate fundamental relation is Sonine's first integral for Bessel functions of the first kind (see [27]). In a recent paper [4] the authors applied similar operators to derive parametrized families of suitable locally supported covariance models for stationary random fields on the sphere S d−1 . These operators are based on properties of Gegenbauer polynomials, appearing in Schoenberg's characterization [24] of zonal positive definite functions on the sphere.
in which all the coefficients a n are nonnegative, and the series converges at the point x = 1. If this is the case, the series converges absolutely and uniformly on the whole interval.
Chen, Menegatto and Sun [7] showed that a necessary and sufficient condition for f • cos to be strictly positive definite on S d−1 , d ≥ 3, is that, in addition to the conditions of Theorem 1.2, infinitely many of the Gegenbauer coefficients a n with odd index, and infinitely many of those with even index, are positive. In the case d = 2 the criteria is necessary but not sufficient for f • cos to be strictly positive definite. A characterization in this case has been given by Menegatto, Oliveira & Peron [18], although the criterion is a little more involved (see also [3] for further details on these issues).
In the same spirit as for the turning bands method, a zonal function defined on a lower dimensional sphere S d−κ can be lifted up to S d−1 through averaging over the set of all copies of S d−κ contained in S d−1 . In [4] it is shown that the analogues for the sphere of Matheron's montée and descente operators (see [16]) for R d are the operators and Paralleling the behaviour of Matheron's operators in the Euclidean case the operators move in the ladder of dimensions by steps of two. Specifically, the I and D operators map zonal positive definite functions f (cos ·) on S d onto ones on S d−2 and S d+2 , respectively (see [4] for details). Therefore, the natural question arises, whether it would also be possible to proceed through steps by one within the ladder of dimensions. While in the Euclidean case this could be achieved using fractional differentiation and integration (see [27]), the situation is more intriguing in the spherical setting. The reason lies in the fact that the characterizing special functions for the sphere are polynomials, which are not preserved through fractional integration. Thus, one has to work with combinations of fractional operators in order to guarantee that the operators are mapping into the space of polynomials.
In the present paper, we provide a suite of four operators which can be used to define a clavier (cf. [16]) for the sphere. The main results are given in Theorems 2.3, 2.4 and 2.8, below. We start with introducing the appropriate fractional operators in the following section and studying their action on ultraspherical expansions. The action of the operators on Gegenbauer polynomials shown in the last section is derived using properties of hypergeometric 2 F 1 -functions.

Def inition of the half-step operators
In the expansion (1.1) the dimension d appears in the parameter λ = (d−2)/2 of the Gegenbauer polynomials. This relation between λ and d will be fixed throughout the paper.
From D C λ n = 2λC λ+1 n−1 (cf. [9, 10.9(23)]) we see that classical differentiation and its inverse, integration, alter the parameter λ by an integer. This is why the operators I and D traverse the ladder of dimensions in steps of two (see [4]). At the same time, I and D change the degree of polynomials by one. Therefore, in order to obtain a one-step operator in the dimension, we have to consider fractional integration and differentiation, a fact which perfectly parallels the Euclidean setting (see [20,27]).
We are now ready to define the half-step operators and discuss their action on positive definite functions on S d−1 .
Using these, we further define Apart from the additional factor (1 ± x) −λ+ 1 2 in front of the integral and the weight (1 ± τ ) λ , the operators I λ ± are classical Riemann-Liouville fractional integrals of order Using these, we further define The main results of this paper are the following two theorems giving precise statements of the dimension hopping and positive definiteness preserving properties of the operators I λ ± and D λ ± . These are one-step analogues of Theorems 2.2 and 2.3 in [4]. Since in the light of Theorem 1.2 the statements can be considered as statements concerning ultraspherical expansions without referring back to a sphere, we are considering m = d − 1 to be a positive integer.
with all coefficients, {a n } ∞ n=0 , positive. Then I λ + f is also nonnegative, I λ + f ∈ Λ + m , and all the coefficients b n in the expansion are positive.
Proof . The proofs for the statements are almost identical with those of the corresponding parts of Proposition 2.2 in [4], provided that proper analogues for certain statements on Gegenbauer polynomials are given. We therefore restrict ourselves to pointing out where adaptations of the proof given in [4] are needed. One of these details concerns the boundedness of the operators I λ ± as operators from C[−1, 1] to C[−1, 1]. This follows from the definitions of I λ ± and I λ ± in equations (2.1), (2.2) and (2.3), combined with the beta integrals , .
Similarly, positivity of the operator I λ + follows from the definitions (2.1), (2.2) and (2.3). The main ingredient thus remaining to be shown is the action of I λ ± on the Gegenbauer polynomial C λ+ 1 2 n . This part is given in Theorem 3.3, below. Note that in contrast to the operators studied in [4], there is no need to deal with an extra constant in statements (i) and (ii). This follows from Theorem 3.3, showing that the operators I λ ± do not introduce an additional constant.
Since the operators defined above can be seen as standard operators of fractional integration/differentiation, classical results from fractional calculus can be applied. For example, if (1 + τ ) −1/2 f (τ ) ∈ Lip α for some α > 1 2 , in particular, if f ∈ Lip α and suppf ⊂ (−1, 1], then by Theorem 19 in [14] D 0 ± f exists and is continuous. An analogous statement holds for general λ. The proof of Theorem 2.4 depends heavily on a multiplier relationships between the Gegenbauer coefficients of f and those of D λ ± f . The details of these relationship, and the proof of Theorem 2.4, will be deferred to the next subsection.
Let us finish the section with considering an example.
Example 2.6. Consider the operator I λ + . In view of its definition (2.1) this operator maps functions locally supported near one to functions locally supported near one. Also, since I λ + = (I λ + + I λ − )/2 this operator preserves (strict) positive definiteness by Theorem 2.3. Note that by a change of variables Therefore, if f were such that f (x + 1)s − 1 = (1 − ys) −a , the integral becomes which is strictly positive definite on R d for all d ≥ 1 (cf. [13]). Choosing α = 2 and restricting the function ϕ 2,β to the sphere, we obtain (setting τ = cos θ) Thus, Therefore, from (2.5) we have that we can choose β = 2λ + 3, yielding Therefore, the function given in (2.7) is strictly positive definite on S 2λ+1 by the remark at the start of the example. We can follow the same line of argument applying I µ + to I λ + ϕ 2λ+3 . Again, the result is a hypergeometric function. Interestingly, for µ = λ − 3 2 the series is of the form (2.6) resulting in In general, for a function

Ultraspherical expansions of f and D λ ± f
The main results of this section will be Theorems 2.7 and 2.8 giving multiplier relationships between the Gegenbauer coefficients of the (formal) series of f and those of the (formal) series of D λ + f and D λ − f . These relationships will later be used to show that the operators D λ + and D λ − have the positive definiteness preserving properties given in Theorem 2.4. The first statement shows that the operators D 0 + , D 0 − can be applied term by term to a Chebyshev series to obtain the formal Legendre series of D 0 ± f . Similar relations between coefficients in ultraspherical expansions hold for higher order Gegenbauer polynomials. If both functions D λ ± f ∈ C[−1, 1], then the (formal) Gegenbauer series 8) and the (formal) Gegenbauer series has coefficients n + λ a n , n ∈ N 0 . (2.9) Remark 2.9. Theorem 2.7 is the limiting case of Theorem 2.8 under the limit Furthermore, we have the special cases Before proving the theorems we state the following technical lemma.
Proof . Formula (2.12) can be obtained from equation (2.11) by using the reflection formula for Gegenbauer polynomials . (2.13) and the change of variables y = −x.
For the proof of formula (2.11) we will use two recurrences involving derivatives of Gegenbauer polynomials which can be found, for example in [9, 10.9(35)]. For notational convenience we use the (non-standard) notation D λ n (x) = d dx C λ n (x) within the proof of the lemma. Then, nC λ n (x) = xD λ n (x) − D λ n−1 (x), and (2.14) Turn now to an induction proof of formula (2.11). The statement is clearly true for C λ 0 (x) = 1, adopting the convention that the sum then is empty. Now assume that n ∈ N and that the first statement is true for n − 1. Consider Using (2.15) and then (2.14) we obtain that Therefore, applying (2.16), Using the induction hypothesis gives which completes the proof.
Proof of Theorem 2.7. Note that the continuity of the two functions D 0 ± f implies that of the functions D 0 ± f . Then, proceeding by integration by parts, the coefficient b n of P n in the formal Legendre expansion of D 0 + f is (2n + 1)/2 times In view of the formula [19, (5.12.1)], for all −1 < x ≤ 1, the limit as x tends to −1 of the quantity within the square brackets vanishes. Hence, where the last step follows from an application of formula (2.11 after some straightforward calculation, the double integral above turns into the form Therefore, Analogously, we define A similar integration by parts argument, but now using the formula (2.12), and the relationship (cf. [19, (18.17.45 shows that H − = − n 2n + 1 πa n + n + 1 2n + 1 πa n+1 . Since b n = 2n+1 2 (H + + H − ) we finally obtain b n = (n + 1)πa n+1 , n ∈ N 0 .
This completes the proof of the part of the theorem concerning D 0 + f . The proof of the part of the theorem concerning D 0 − f is similar and will be omitted.
The proof of Theorem 2.8 relies on a kind of fractional integration by parts. Before going into details, we will state some technical lemmas.
The following proposition is the limit case of (2.20) taking the limit λ → 1 + after multiplying either side with 1/(λ − 1).
Proof . From [9, 10.11(3)] we have that Furthermore, [9, 10.11(37)] yields The claim follows from using these relations to rewrite the right hand side in terms of Chebyshev polynomials of the first kind.
The following lemma states the analogues of equations (2.17) and (2.18) for general λ > 0. Note that the first integral given in the lemma is a special case of Bateman's integral (see [19, (18.17.9)]). Using this, the statement follows from [19, (18.9.4)]. Rather than using Bateman's integral for Jacobi polynomials, we provide a proof staying within the family of Gegenbauer polynomials, only.
Lemma 2.14. Let λ > 0, n ∈ N 0 and x ∈ [−1, 1]. Then Proof . From the definition we have that and 2I λ − = I λ + − I λ − . We can therefore use Theorem 3.3 below to compute from which the first statement follows. The second integral follows from the first by a change of variables τ = −t and using the reflection formula (2.13).
Proof of Theorem 2.8. The Fourier-Gegenbauer coefficients of the function D λ + f are given by Let us denote the first integral by H + and the second by H − . From integration by parts it follows that The constant term vanishes. Applying Lemma 2.11 we can decompose H + into the sum of the two integrals and In an analogous way, we decompose H − into a sum of the integrals and The inner integrals in I 1 and I 3 can be computed using Lemma 2.14, whereas the corresponding formulae for I 2 and I 4 are stated in Corollary 2.16. Doing so, and using the definition of the coefficients we obtain that 2 2λ n!Γ(λ)Γ λ + 1 2 n + λ + 1 2 × Γ(2λ + n + 1) λ + n + 1 a n+1 − (2λ + n) Γ(2λ + n) λ + n a n , × (2λ + n) Γ(2λ + n + 1) λ + n + 1 a n+1 + (2λ + n) Γ(2λ + n) λ + n a n , × (n + 2λ) Γ(2λ + n) λ + n a n + Γ(2λ + n + 1) λ + n + 1 a n+1 , Therefore, from which it follows that This completes the proof of the part of the theorem concerning D λ + f . Again, the proof of the part concerning D λ − f is similar and will therefore be omitted.
Proof of Theorem 2.4. From the continuity assumption we have that the Gegenbauer series of the functions D λ ± f are Abel summable. Since f is positive definite by hypothesis, Theorem 2.8 shows that the Gegenbauer coefficients of D λ ± f are non-negative. Hence the Gegenbauer series restricted to x = 1 are series of non-negative terms. For such series of non-negative terms Abel summability implies summability. Since C 3 The action of I λ ± and D λ ± on Gegenbauer polynomials The claim that the operators I λ ± and D λ ± map Gegenbauer polynomials onto Gegenbauer polynomials with a changed parameter is based upon contiguous relations for hypergeometric functions.
We will thus first state a proposition showing that the images of C   (1 + x) 2 F 1 −n, n + 2λ + 1 λ + 1 .
Proof . Consider the first equation Substituting into equation (3.12) gives equation (3.10). The proof of equation (3.11) is almost identical to that of equation (3.10). It will therefore be omitted.
The proof of equation (3.9) is analogous except that it uses contiguous relation [19, (15.5.16)]. Hence, it will be omitted.