Abstract

A new reweighted proportionate affine projection algorithm (RPAPA) with memory and row action projection (MRAP) is proposed in this paper. The reweighted PAPA is derived from a family of sparseness measures, which demonstrate performance similar to mu-law and the l0 norm PAPA but with lower computational complexity. The sparseness of the channel is taken into account to improve the performance for dispersive system identification. Meanwhile, the memory of the filter’s coefficients is combined with row action projections (RAP) to significantly reduce computational complexity. Simulation results demonstrate that the proposed RPAPA MRAP algorithm outperforms both the affine projection algorithm (APA) and PAPA, and has performance similar to l0 PAPA and mu-law PAPA, in terms of convergence speed and tracking ability. Meanwhile, the proposed RPAPA MRAP has much lower computational complexity than PAPA, mu-law PAPA, and l0 PAPA, etc., which makes it very appealing for real-time implementation.

Keywords

1 Introduction

Adaptive filtering has been studied for decades and has found wide areas of application. The most common adaptive filter is the normalized least mean square (NLMS) algorithm due to its simplicity and robustness [1]. In the 1990’s, the affine projection algorithm (APA), a generalization of NLMS was found to have better convergence than NLMS for colored input [2, 3]. The optimal step size control of the adaptive algorithm has been widely studied in order to improve their performance [4, 5]. The impulse responses in many applications, such as network echo cancellation (NEC), are sparse, that is, a small percentage of the impulse response components have a significant magnitude while the rest are zero or small. To exploit this property, the family of proportionate algorithms was proposed to improve performance in such applications [2]. These algorithms include proportionate NLMS (PNLMS) [6, 7], and proportionate APA (PAPA) [8], etc.

The idea behind proportionate algorithms is to update each coefficient of the filter independently of the others by adjusting the adaptation step size in proportion to the magnitude of the estimated filter coefficient [6]. In comparison to NLMS and APA, PNLMS and PAPA have very fast initial convergence and tracking when the echo path is sparse. However, the big coefficients converge very quickly (in the initial period) at the cost of slowing down dramatically the convergence of the small coefficients (after the initial period). In order to combat this issue, mu-law PNLMS (MPNLMS) and mu-law PAPA algorithms were proposed [9–11]. Furthermore, the l0 norm family of algorithms have recently drawn lots of attention for sparse system identification [12]. Therefore, a new PNLMS algorithm based on the l0 norm was proposed to represent a better measure of sparseness than the l1 norm in PNLMS [13].

On the other hand, the PNLMS and PAPA algorithms converge much slower than corresponding NLMS and APA algorithms when the impulse response is dispersive. In response, the improved PNLMS (IPNLMS) and improved PAPA (IPAPA) were proposed by introducing a controlled mixture of proportionate and non-proportionate adaptation [14, 15]. The IPNLMS and IPAPA algorithms perform very well for both sparse and non-sparse systems. Also, recently, the block-sparse PNLMS (BS-PNLMS) algorithm was proposed to improve the performance of PNLMS for identifying block-sparse systems [16].

In order to reduce the computational complexity of PAPA, the memory improved PAPA (MIPAPA) algorithm was proposed to not only speed up the convergence rate but also reduce computational complexity by taking into account the memory of the proportionate coefficients [17]. Dichotomous coordinate descent (DCD) iterations have previous been applied to the PAPA family of algorithms to implement the MIPAPA adaptive filter [18, 19]. Meanwhile, an iterative method based on the PAPA with row action projection (RAP) has been shown to have good convergence properties with relatively low complexity [20].

In [21] the proportionate adaptive filter was derived from a unified view of variable-metric projection algorithms. In addition, the PNLMS algorithm and PAPA can both be deduced from a basis pursuit perspective [22, 23]. A more general framework was further proposed to derive PNLMS adaptive algorithms for sparse system identification, which employed convex optimization [24]. Here, a family of PAPA algorithms are firstly derived based on convex optimization, in which PAPA, mu-law PAPA, and l0 PAPA are all special cases. Then, a reweighted PAPA is suggested in order to reduce the computational complexity. Finally, an efficient implementation of PAPA is proposed based on RAP and memory PAPA.

The organization of this article is as follows. The review of various PAPAs is presented in Section 2. Section 3 derives the proposed reweighted PAPA and presents an efficient memory implementation with RAP. The computational complexity is compared with PAPA, mu-law PAPA and l0 PAPA in Section 4. In Section 5, simulation results of the proposed algorithm are presented. The last section concludes the paper with remarks.

2 Review of various PAPAs

The input signal x(n) is filtered through the unknown coefficients to be identified h(n) to get the observed output signal d(n).

where \(\mathrm {F}(|\hat {h}_{l}|)\) is specific to the algorithm, q prevents the filter coefficients \(\hat {h}_{l}(n-1)\) from stalling when \(\hat {\mathbf {h}}(0)=\mathbf {0}_{L\times 1}\) at initialization and ρ prevents the coefficients from stalling when they are much smaller than the largest coefficient. The classical PAPA employs step-sizes that are proportional to the magnitude of the estimated impulse response as below [8]

Based on the motivation that the l0 norm can represent an even better measure of sparseness than the l1 norm, the improved PNLMS and PAPA algorithms based on an approximation of the l0 norm (l0-PNLMS) were proposed as below [13]:

where σl0 is a positive parameter. The main disadvantage of the mu-law or l0 norm PAPA algorithms are their heavy computation cost because of the L logarithmic or exponential operations. Therefore, a line segment was given to approximate the mu-law function [9], where

It should be noted that, without loss of performance, the line segment was normalized to be of unit gain for \(|\hat {h}_{l}|\geq 0.005\), compared to the original one proposed in [9]. Meanwhile, the exponential form in (12) can be approximated by the first order Taylor series expansions of exponential functions [12]

It is interesting to see that the first order Taylor series approximation of l0 PAPA in (12) is actually the same as the line segment implementation of mu-law PAPA in (11) for σl0=200.

3 The proposed SC-RPAPA with MRAP

Based on the minimization of the convex target, the reweighted PAPA (RPAPA) will be firstly derived from a new sparseness measure with low computational complexity. Meanwhile, the sparseness controlled RPAPA (SC-RPAPA) is presented to improve the performance for both sparse and dispersive system identification. Finally, the SC-RPAPA with memory and RAP (MRAP) is proposed by combing the memory of the coefficients with iterative RAP to further reduce the computational complexity.

3.1 The proposed RPAPA

The proportionate APA algorithm can be deduced from a basis pursuit perspective [22]

where G−1(n) is the inverse matrix of proportionate matrix G(n), which is also a diagonal matrix. If the optimization target in (17) is convex, the family of PAPA algorithms can be derived using Lagrange Multipliers. It should be noted that, using the approximation

the proposed formulation in (17) becomes the variable-metric in [21], which is an approximation of the proposed formulation. The function \(\mathrm {G}(t), t\in \mathbb {R}\) should satisfy the following properties:

1) G(0)=0, G(t) is even and not identically zero;

2) G(t) is non-decreasing on [0,∞);

3) \(\frac {\mathrm {G}(t)}{t}\) is non-increasing on (0,∞).

The above properties follow the requirements of the sparseness measure proposed in [25]. From the perspective of proportionate algorithms, the first two requirements are intuitive, since the family of the proportionate algorithms should be proportionate to the magnitude of the filter’s coefficients. The third property will guarantee the convexity of the optimization target. PAPA, mu-law PAPA and l0 PAPA are all special cases of the sparseness measures fulfilling all three properties. In this paper, considering the computational complexity, we propose using the following reweighted PAPA:

The proposed reweighted metric is compared with PAPA, mu-law PAPA and l0 PAPA in Fig. 1. The σ parameters for each algorithm were σμ=1000, σl0=50, σr=0.01. These parameters were recommended and widely simulated in the literature for each algorithm [9, 13]. It should be noted that, the plots in [24] set the σ parameters respectively so that they all contain the point (0.9,0.9). However, in actual application, this parameter should be tuned to maximize the performance. In order to facilitate the comparison of the different sparseness measure, they are normalized to pass through the point, (1,1) here instead. Without loss of generality, it is assumed that the filter’s coefficients are normalized and the maximum possible magnitude is 1. Therefore, it is convenient to compare the gain distribution of different metrics with different σ parameters.

Fig. 1

Comparison of the different metrics

3.2 The proposed SC-RPAPA

It should be noted that the reweighting factor σr in the proposed RPAPA (19) is related to the sparseness of the impulse system. It is straightforward to verify that if σr=0, reweighted PAPA simplifies to APA. If the impulse system is more sparse, σr should be relatively larger than \(\left |\hat {h}_{l}\right |\), which makes it more like the PAPA. This agrees with the fact that we fully benefit from PNLMS only when the impulse response is close to a delta function [26]. Therefore, it is natural to take the sparseness of impulse response into account. The sparsity of an impulse response could be estimated as

where L>1 is the length of the channel, \(\Vert \hat {\mathbf {h}}(n)\Vert _{1}\) and \(\Vert \hat {\mathbf {h}}(n)\Vert _{2}\) are the l1 norm and l2 norm of \(\hat {\mathbf {h}}\), respectively. The value of \(\hat {\epsilon }(n)\) is between 0 and 1. For a sparse channel, the value of the sparseness is close to 1 and for a dispersive channel, this value is close to 0. Therefore, the SC-RPAPA is

where σmax is the maximum value for the sparse system identification. The plot of the reweighted metric for different σs is presented in Fig. 2. In practical implementation, we would like to apply the APA algorithm to the dispersive system under certain sparseness threshold. For example, the sparsity of the dispersive channel is about 0.4, and a heuristic implementation that works pretty well in the simulations is

where εmin=1e−4 is a minimum sparsity in order to avoid dividing by zero for \(\hat {h}_{l}=0\).

3.3 The proposed SC-RPAPA with MRAP

However, the main computational complexity of the family of PAPA algorithm is the matrix inversion in (5). Reduction in complexity is achieved by using 5M DCD iterations, thus requiring about 10M2 additions [18]. Meanwhile, a sliding-window recursive least squares (SRLS) low-cost implementation of PAPA is given based on DCD, which does not depend on M. The SRLS implementation is only efficient when the projection order is very high (e.g., such as M=512) [19]. However, it is known that if the projection order increases, the convergence speed is faster, but the steady-state error also increases.

Another way to avoid the matrix inversion altogether is to use the method of RAP [27]. RAP is also known in the literature as a data reuse algorithm (see [28]). It has been shown in [29] that RAP is effectively the same as APA, except that the system of equations problem that is solved with a direct matrix inversion (DMI) in APA is solved iteratively in RAP [20].The iterative PAPA algorithm proposed in [30] was made efficient by implementing it using RAP in [27]. RAP is an iterative approach to solving a system of M equations. It cycles through the M equations J times performing an NLMS-like update on the coefficients for each equation. In this instance, the number of RAP iterations, J is set to one. It should be noted that, by limiting J to one, the solution of the system of equations through RAP is approximate. However, the simulation results will demonstrate that this approximation works pretty well, especially for relatively high projection order. In each sample period a new equation is added to the system of equations and the oldest equation is dropped. Thus, M RAP updates are performed on a given equation every M sample periods. The PAPA algorithm with RAP updates the coefficients

in which the operation ⊙ denotes the Hadamard product and m=0,1,…,M−1.

The traditional PAPA requires M×L multiplications to calculate P(n), and in order to further reduce the computational complexity, we propose to apply the memory of the proportionate coefficients [17] into SC-RPAPA. Therefore, the matrix P(n) in (4) can be approximated as P′(n)

As mentioned in [17], the proposed RPAPA with MRAP takes into account the “history” of the proportionate factors from the last M steps. The convergence and the tracking become faster when the projection order increases. Meanwhile, combined with the RAP, the computational complexity is also significantly lower as compared to the MPAPA through avoiding the direct matrix inversion and using the memory. The proposed SC-RPAPA with MRAP algorithm is summarized in detail in Table 1.

4 Computational complexity

The computational complexity of the SC-RPAPA with MRAP algorithm is compared with traditional PAPA, MPAPA, RPAPA, and SC-RPAPA in Table 2, in terms of the total number of additions (A), multiplications (M), divisions (D), comparisons (C), square root (Sqrt), and direct matrix inversion (DMI) needed per algorithm iteration. All the algorithms require L |·| operations for calculating the magnitude of the filter’s coefficients.

Table 2

Computational complexity of the algorithms’ coefficient updates

A

M

D

C

Sqrt

DMI

PAPA

(M2+2M+1)L−M−1

(M2+3M+1)L+2M2+2

L

2L

0

Yes, M×M

MPAPA

(M2+2M+1)L−M−1

(M2+2M+2)L+2M2+2

L

2L

0

Yes, M×M

RPAPA

(M2+2M+2)L−M−1

(M2+3M+1)L+2M2+2

2L

2L

0

Yes, M×M

SC-RPAPA

(M2+2M+4)L−M−1

(M2+3M+2)L+2M2+5

2L+1

2L+1

1

Yes, M×M

SC-RPAPA MRAP

(2M+5)L+M−2

(2M+3)L+M+5

2L+M+1

2L+1

1

No

Compared with traditional PAPA, the MPAPA reduced the complexity of GX, but the calculation of XTP′ still requires M2L multiplications. Meanwhile, due to the memory and the iterative RAP structure, only L multiplications are needed to update p(n) instead.

What’s more important is that, both the PAPA and the MPAPA algorithms require a M×M direct matrix inversion, which is especially expensive for high projection orders. The combination of the memory and the iterative RAP structure, not only avoids the M×M direct matrix inversion, but also reduces the computational complexity required for the calculation of both GX and XTGX.

The additional computational complexity for the SC-RPAPA with MRAP algorithm arises from the computation of the sparseness measure \(\hat {\epsilon }\). As in [31], given that \(L/(L-\sqrt {L})\) can be computed offline, the remaining l-norms require an additional 2L additions and L multiplications. Furthermore, this sparseness measure can be reused in many other sparseness controlled algorithms too, for example [31]. The calculation of the F in (22) requires additional L divisions, L+1 additions, one multiplication, and one comparison more than PAPA. The complexity of division is much lower than the L exponential or logarithmic operations required by either the mu-law or the l0 PAPA. Meanwhile, (22) also offers the robustness to dispersive system identification.

5 Simulation results

The performance of the proposed SC-RPAPA with MRAP was evaluated via simulations. Throughout our simulation, the length of the unknown system was L=512, and the adaptive filter was with the same length. The sampling rate was 8 kHz. The parameters for each algorithm were δ=0.01/L, ρ=0.01, q=0.01. The step-size for all the algorithms was set to μ=0.2.

The algorithms were tested using both the white Gaussian noise (WGN), and colored noise as inputs. The colored input signals were generated by filtering the WGN through a first order system with a pole at 0.8. Independent WGN was added to the system background with a signal-to-noise ratio (SNR) as 30 dB.

Two impulse responses were used to verify the performance of the proposed SC-RPAPA MRAP algorithm, as shown in Fig. 3. The first one in Fig. 3a is a sparse impulse response of typical network echo with sparseness 0.92. Figure 3b is a dispersive channel with sparseness 0.44. In order to demonstrate the tracking ability, an echo path change was incurred through switching the impulse response from the sparse system in Fig. 3a to the dispersive one in Fig. 3b.

Fig. 3

Two impulse responses used in the simulation a the sparse network echo path, and b the dispersive echo path

The convergence state of adaptive filter is evaluated with the normalized misalignment which is defined as

5.1 The performance of the proposed RPAPA

The proposed reweighted PAPA in (19) was firstly compared to PAPA, mu-law PAPA, and l0 PAPA. The parameters for the algorithm were σμ=1000, σl0=200, and σr=0.01. The affine projection order was selected as M=2.

In the first simulation shown in Fig. 4, the input signal was the WGN. According to the results, the proposed RPAPA could outperform PAPA, and has similar performance with respect to mu-law and l0 PAPA. However, the reweighted PAPA has much lower computational complexity. In the second simulation, the input signal was colored, and a similar result could be obtained according to Fig. 5.

5.2 The performance of the proposed SC-RPAPA

To demonstrate the benefit of sparseness control, the proposed SC-RPAPA algorithm was simulated using an echo path change from the sparse to the dispersive impulse response in Fig. 3. The SC-RPAPA algorithm was compared with APA, PAPA, and the above RPAPA algorithms. The parameters for the algorithm were σr=0.01, and σmax=0.02. The affine projection order was selected as M=2. In Fig. 6, the input signal was the WGN input. Both the proposed RPAPA and SC-RPAPA algorithms had similar performance for sparse system identification, which outperformed APA and PAPA. Meanwhile, due to the sparseness control, SC-RPAPA outperformed RPAPA as expected for the dispersive system. The colored input was used in Fig. 7, and similar results are observed.

5.3 The performance of the proposed SC-RPAPA with MRAP

An efficient implementation of the SC-RPAPA algorithm was proposed through combining the memory of the filter’s coefficients with RAP. The new SC-RAPA with MRAP algorithm significantly decreases computational complexity. In this subsection, the performance of the efficient implementation was compared with APA, PAPA and SC-RPAPA through simulations.

In the first simulation, the WGN input was used. As shown in Fig. 8, SC-RPAPA with MRAP worked as well as SC-RPAPA for sparse system identification. However, for dispersive system, the performance of SC-RPAPA MRAP was worse than SC-RPAPA and the APA. This fact becomes more apparent for the colored input as shown in Fig. 9. This was caused by the relatively low projection order (M=2), and the implementation of the MRAP was slower than the direct matrix inversion. However, this drawback could be mitigated through increasing the projection order. Furthermore, the memory of the filter’s coefficients will also improve the performance as the projection order increases. We verify this point through simulations with M=32 for both the WGN (see Fig. 10) and the colored input (see Fig. 11). It could be observed that the SC-RPAPA with MRAP works better than APA, PAPA, and SC-RPAPA for sparse system identification. Meanwhile, the performance for dispersive system with colored input has been significantly improved too.

6 Conclusion

A low complexity reweighted proportionate affine projection algorithm was proposed in this paper. The sparseness of the channel was taken into account to improve the performance for dispersive systems. In order to reduce computational complexity, the direct matrix inversion of PAPA was iteratively implemented with RAP. Meanwhile, the memory of the filter’s coefficients were exploited to improve the performance and further reduce the complexity for high projection orders. Simulation results demonstrate that the proposed sparseness controlled reweighted proportionate affine projection algorithm with memory and RAP outperforms traditional PAPA, with much lower computational complexity compared to mu-law and l0 PAPA.

Declarations

Acknowledgements

This work was performed under the Wilkens Missouri Endowment. The authors would like to thank the Associate Editor and the reviewers for their valuable comments and suggestions.

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.

Competing interests

The authors declare that they have no competing interests.

Authors’ information

Steven L. Grant was formerly Steven L. Gay.

Authors’ Affiliations

(1)

Department of Electrical and Computer Engineering, Missouri University of Science and Technology, Rolla, USA