Causal Restriction and Its Generalization for the Wiener Filter

The classical principles of Wiener filter design consider only the causality constraint. However, additional constraints are often involved in practice. This paper provides a comprehensive theoretical approach for designing Wiener filters with constraint conditions beyond causality. Using integration, we first derive an extension of the Wiener–Hopf equation with an additional term that captures the effect of the constraint conditions on the resulting Wiener filters. Next, we analytically solve the extended Wiener–Hopf equation and develop a new iterative algorithm to obtain its optimum value. We prove the convergence of our iterative algorithm and provide several simulations to illustrate the applicability of our theoretical approach.


Introduction
The Wiener filter [15] is a well-known tool in signal processing and communications [9]. It presents a classical approach to determining the optimum filter by providing  the minimum mean-squared error estimate of a desired process d (t), from a noisy random input signal x (t), as illustrated in Fig. 1.
In Fig. 1, x (t) and d (t) are assumed to be jointly wide-sense stationary random processes with known autocorrelations, R x (τ ) = E {x(t + τ )x * (t)} and R d (τ ) = E {d(t + τ )d * (t)}, and a known cross-correlation, R dx (τ ) = E {d(t + τ )x * (t)} , in which E {•} denotes the mathematical expectation and the asterisk indicates the complex conjugation operation. The output of the filter isd(t), which denotes the estimate of d (t) and is given bŷ The estimated error is The mean-squared error is The optimum filter, denoted by w opt (t), is determined by minimizing E |e(t)| 2 .
In the classical principle of the Wiener filter design, there is only one causality restriction on the filter to be considered. However, in practice, to meet the needs of different applications, the filter may be subject to more broad constraints. For example, in signal prediction problems and adaptive control applications, the filter is generally required to be a time-limited system to realize finite-step prediction or finite-order control [10]. In past years, numerous researchers have studied the design of the Wiener filter with finite support in the time domain [4,6,14]. Other common restrictions on the Wiener filter are in the frequency domain. For example, the Wiener filter may be required to be band-limited, band-pass, or low-pass or to satisfy other bandwidth limitations to guarantee that the frequency support of the filter can match the bandwidth of the desired signal and to avoid situations in which the filter works in noise-only periods [7,8,12]. For instance, in a speech enhancement system, to reduce the effects of environmental noise, the Wiener filter usually operates in the audible frequency range. Much research has been reported on constrained Wiener filters in binaural hearing aids and speech enhancement [1,2,4,5,12]. In addition to time-domain and frequencydomain constraints, linear constraints, which require the Wiener filter to satisfy a set of linear equations, are also widely applied in many fields, including adaptive beamforming [13], spectral analysis, and function interpolation [11]. Moreover, in filter design, when the power spectrum of the input is close to zero at some frequency points, the corresponding frequency response of the filter at these frequency points will become unbounded. A filter with a boundless frequency response is undesirable because it will be sensitive to noise and data errors. In this situation, the finite-energy constraint is imposed to smooth the peak value of the frequency response of the filter. The performance of the Wiener filter with the finite-energy constraint is better than that of the unconstrained Wiener filter [3]. In addition to the common constraints described above, studies on wider classes of constraints for the Wiener filter can be found in the literature [7,8], including the output-spectrum constraint, the estimated signal-power constraint, the filter-energy constraint, and the filter-stability constraint. Moreover, the closed-form solutions to these constrained design problems have been derived.
The work presented in this paper makes three contributions. First, this is the first work to incorporate the different constraint conditions into a single integration equation. Moreover, since the functions in the integral can be chosen according to different application needs, single integration allows more flexibility and applicability in processing various restrictions on the Wiener filter. Second, we develop a closed formula to represent the optimal solution of the Wiener filter under the constrained integration equation. Because the closed formula has one additional term compared to the Wiener-Hopf equation, the formula is called the extended Wiener-Hopf equation in this paper. Third, we develop a uniform analytical method along with an iterative algorithm to solve the extended Wiener-Hopf equation. Finally, some simulation examples are given to further illustrate the solution steps of our algorithm.

Constraint Condition
Assume that w(t) is subject to the following constraint: where g(t), λ(t), and β may be fully known or some of their characteristics may be known.
In particular, Eq. (4) is a more general expression of the following three constraints.

Time-Domain Constraints
In this case, w(t) is equal to a given function g(t) within the time subset T , i.e., More specifically, if g(t) = 0 and T = (−∞, 0), then w(t) is obviously required to be causal; if g(t) = 0 and T = (−∞, 0) ∪ [t 0 , ∞), where t 0 > 0, then w(t) is a causal time-limited system.

Orthogonal Constraints
In this case, w(t) is assumed to be orthogonal to the signal space H, i.e., Equation (8) is a special case of Eq. (4), in which g(t) = 0 and β = 0. In practice, the orthogonal space of the filter w(t) is usually viewed as noisy space or undesirable signal space, which needs to be depressed by introducing the constraint defined in Eq. (8).
According to these three particular cases, we can construct the Wiener filter with more general constraints to meet the requirements of different practical applications because there are many possibilities for g(t) or λ(t).

Principle
To consider the constraints of Eq. (4), we first multiply the two sides of Eq. (4) by a Lagrange factor μ. As a result, we have μ * ( Then, we add it and its conjugate to Eq. (3) as follows: , and λ(t), respectively. Note that the Fourier transformation of , then L (w, μ) achieves the minimum of E |e(t)| 2 , denoted by L min w opt , μ , and we have and or, equivalently, Correspondingly, by the inverse Fourier transform, in the time domain, Eq. (13) can be converted into a convolution integral of the following form, Equation (14) has one more term, −μλ (τ ), than the Wiener-Hopf equation [9], i.e., , which is easily derived from the unconstrained Wiener filter design problem [9] or from Eq. (9) or (10) by ignoring the integral of the constraint. Since μλ (τ ) is the function of the constraint on the design result and since μλ (τ ) = 0 if and only if there is no constraint on the Wiener filter, Eq. (14) can be called the extended Wiener-Hopf equation. Finally, we take our familiar causal constraint as an example to further illustrate the role of μλ (τ ) in the extended Wiener-Hopf equation. As stated in the discussion in Eq. (5), requiring the filter to be causal is equivalent to requiring that (4). Consequently, in the corresponding extended Wiener-Hopf equation under the causal constraint, μλ (τ ) is an anti-causal function whose purpose is to offset the non-causal component of the Wiener-Hopf equation and to guarantee that the optimal filter is a causal system.
In addition, we develop the other computational formulas of the minimum meansquared error of the Wiener filter under the integral constraint. From Eqs. (3), (12), (13), and (9), we can obtain The first equality is immediately derived from the Fourier transformation of Eq. (3), after w(t) is replaced by w opt (t). The second equality is obtained by replacing The third equality is obtained by replacing S in the second equality. The last equality is a simplification of the third equality since μ

Algorithms
In this section, we discuss how W opt ( f ) can be obtained from Eq. (12). First, we give the solution of Eq. (12) for the constraint condition in which the filter is required to have finite frequency support. Then, we develop an analytical approach and an iterative technique to obtain W opt ( f ) from Eq. (13).
In practice, the filter may be required to be low-pass, band-pass, high-pass, or band-limited. These restrictions can be summarized as the finite frequency-support constraint, which is explained as follows: If the filter is required to satisfy where denotes a non-empty subset of the frequency domain, then the filter has finite frequency support. In this case, ( f ) should be chosen such that Then ,we have the corresponding constrained integral From Eqs. (17) and (18) Then, we multiply both sides of Eq. (12) by Γ ( f ). As a result, we obtain the solution of the Wiener filter with finite frequency support: Next, we introduce the analytical approach. By inserting Then, solving this integration equation, we obtain the Lagrange factor μ as As a result, W opt ( f ) can be analytically calculated from Eqs. (20) and (12) if the analytical form of ( f ) is known. However, this approach is invalid when ( f ) is unknown or when only its characteristics are known without any information about its form. Now, we introduce an iterative technique, which is valid regardless of whether the analytical form of ( f ) is known. The iterative formulas are given in the frequency domain by Eq. (21): where η is the convergence factor and the subscript (i) denotes the iteration number. W i+1 ( f ) represents the result of iteration (i+1), and W i+1 ( f ) is the largest constrained component of W i+1 ( f ), which will be defined in the next paragraph. Because W i+1 ( f ) generally does not meet the constraint defined in Eq. (4), for the sake of the restriction condition in the next iterative operation and fast convergence, we need to extract the largest constrained component of W i+1 ( f ) to satisfy Eq. (4) . The iterative procedure can be summarized as follows: (1) Set the convergence factor 0 < η < 2 max S x ( f ) and a positive number ε, which is the iteration termination threshold. Then, arbitrarily choose an initial value W 0 ( f ) and insert it into Eq. (21a) to obtain W 1 ( f ).
(2) Using Eq. (21b) through Eq. (21d), or by immediately applying the constraint condition, extract the largest constrained component of Before discussing the convergence, we define the largest constrained component.

Definition 1 The largest constrained component of
, is defined as follows: W i ( f ) satisfies the above definition and thus is the largest constrained component.

(4) and can easily obtain
However, this orthogonality is possible only if ∂ W i ( f ) is equal to zero. As a consequence, W i ( f ) − W i ( f ) satisfies item (c) of the above definition and is thus the largest constrained component of W i ( f ). This completes the proof of existence.
We can rewrite the calculation formula for the largest component of W i ( f ), i.e., W i ( f ), as follows: where In other words, using the result of the ith iteration of W i ( f ) and Eqs. (22) and (23), we can easily obtain W i ( f ). Then, the next iteration will be implemented.
Finally, it should be mentioned that because can be immediately derived from Eq. (23). If only some characteristics of W ( f ) or ( f ) are known, we can make use of these characteristics to extract W i ( f ) from W i ( f ). In these cases, W i ( f ) is sometimes easily extracted; in Sect. 5, we will present a simulation example for further illustration. Now we discuss the uniqueness. Let us assume that the largest constrained component of W i ( f ) is not unique, namely both W i,0 ( f ) and W i,1 ( f ) satisfy the above definition, along with their residuals W i,0 ( f ) and W i,1 ( f ), respectively. It is easy to show that Eq. (4), which means that W i,1 ( f ) cannot be the largest component if ∂ ( f ) = 0, which contradicts the conjecture. Therefore, ∂ ( f ) is equal to zero. As a result, This completes the proof of the uniqueness.
Proof According to the definition, W i+1 ( f ) is the largest constrained component of W i+1 ( f ) and satisfies the restriction condition. Thus, we have Hence, by subtracting Eq. (26) from (27), we obtain In other words, After finishing the orthogonality analysis, we will prove Eq. (25) in the following discussion.
By subtracting W opt ( f ) from both sides of Eq. (21) and substituting S This completes the proof of convergence for the iterative algorithm. When will converge to zero in the same norm. Hence, we can take d (W i+1 , W i ) < ε as the criterion to terminate the iterative procedure. In addition, we define the approximation error of W opt ( f ) by is positively correlated with ε. Therefore, to obtain a high approximation degree, ε should be as small as possible, but this increases the number of iterations that must be computed. Finally, it should be mentioned that although the iterative technique is developed in the frequency domain, it may be possible to implement this technique in the time domain by making use of the Fourier transformation.

Examples
In this section, three simulation examples are presented to illustrate the applicability of our methods. In all three examples, we take the same noisy random input x(t) and the same desired output d(t), but impose different restrictions. We assume that R x (τ ) =

Case 1
The Wiener filter is required to be a low-pass system with a 40-Hz bandwidth, i.e., Obviously, this is a finite frequency-support constraint. From Eqs. (16), (18), and (19), we obtain The wave form of W opt ( f ) is shown in Fig. 2.

Case 2
The Wiener filter is required to satisfy the following integration equation: It is straightforward to show that this integral is equivalent to Eq. (4) when G( f ) = 0 and The computational procedure of the analytical method is very simple. First, using Eq. (20), we obtain μ = 1.36 × 10 −5 . Then, from Eq. (12), we obtain The wave form of W opt ( f ) is shown in Fig. 3. Next, we use the iterative algorithm to solve this problem again. The initial parameters of this iteration is chosen as 0.02+0.01 4 +(0.02π f ) 2 , ε = 10 −6 , and η = 0.0099 < 2 max S x ( f ) ≈ 0.01. In the above iterative algorithm, after 171 iterations, the criterion to end the iterative procedure is satisfied. As a result, we obtain W 171 ( f ) as shown in Fig. 4.
From Eq. (32), we obtain the approximation error, δ W opt ( f ), W i+1 ( f ) = 1.37%. We can further decrease δ W opt ( f ), W i+1 ( f ) by reducing ε, but the computation time will increase. In this instance, if we take ε = 10 −13 and hold the other initial parameters constant, then, after 162,834 iterations, the iterative procedure is terminated, and we obtain δ W opt ( f ), W 162834 ( f ) = 0.33%.

Case 3
The Wiener filter is required to satisfy the following constraint condition: unknown, t ∈ [0, 500) (36) For this restriction, we take λ(t) = unknown, t / ∈ [0, 500) 0, t ∈ [0, 500) Then we can convert the restriction condition to the integration equation Because λ(t) is unknown, we can use only the iterative algorithm. We implement the iterative procedure in the time domain. We set w 0 (t) = I FT S dx ( f ) S x ( f ) , ε = 10 −25 , and η = 0.0099, where I FT {•} denotes the inverse Fourier transform. In this instance, we cannot use Eq. (23) to extract the largest component of w i (t), because λ(t) is unknown. However, according to the constraint given in Eq. (36), we can immediately obtain w i (t) = g(t) ,t / ∈ [0, 500) w i (t) ,t ∈ [0, 500) For this instance, in the fourth iteration, the criterion for terminating the iterative procedure is satisfied. We obtain w 4 (t) as shown in Fig. 5.

Conclusions
This paper presents a general theoretical method for designing Wiener filters with additional constraint conditions. By introducing an integral expression to represent the constraints, we derive and analytically solve the formula for a more constrained Wiener filter. Compared to the classical design method, our new approach not only incorporates broader constraints but also relaxes the requirement of rational power spectra. The classical method of designing Wiener filters is appropriate only for cases of causality, and the causal solution is obtained only if the input self-power spectrum and the cross-power spectrum of the input and the output are both rational functions. Finally, this paper demonstrates several practical applications of our techniques using simulations.