# The Remez Algorithm

The Remez algorithm is a well-known iterative routine to approximate a function by a polynomial in the minimax norm. But, as Nick Trefethen [1] says about it:

Most of these [implementations] go back many years and in fact, most of them do not solve the general best approximation problem as posed above but variants involving discrete variables or digital filtering. One can find a few other computer programs in circulation, but overall, it seems that there is no widely-used program at present for computing best approximations.

One can compute the minimax solution also by applying least-squares or convex optimization, for example using Matlab and the free CVX toolbox applied to the Runge function on [-1, 1]:

m = 101; n = 11;            % 101 points, polynomial of degree 10
xi = linspace(-1, 1, m);    % equidistant points in [-1, 1]
ri = 1 ./ (1+(5*xi).^2);    % Runge function

tic                         % p is the polynomial of degree (n-1)
cvx_begin                   % minimize the distance in all points
variable p(n);
minimize( max(abs(polyval(p, xi) - ri)) );
cvx_end
toc                         % 0.17 sec for Matlab, CVX and SeDuMi


The approximation with Chebyshev polynomials has a minimax norm of 0.1090 while this approach here reaches a minimum of 0.0654, the same value that gets computed with the Remez algorithm in the Matlab chebfun toolbox.

Is there any advantage in applying the more complicated Remez algorithm if you can compute the minimax solution faster and more accurate with an optimization solver? Are there any reports/articles comparing these two approaches on some difficult problems or test cases?

--
[1] R. Pachon and L. N. Trefethen. BIT Numerical Mathematics (2008) Vol. 46.

The "right" answer strongly depends on what you need your approximant for. Do you really need the best approximation for some error bound? Or just a good approximation? Or just a good approximation in the minmax sense?

Nick Trefethen recently gave a nice example where Remez approximation is a bad idea since it minimizes the maximum error irrespective of the average error over the entire interval, which may not be what you want. Of course, the maximum error may be large, but this is bounded for smooth functions.

Update

Following the discussion in the comments below, I downloaded the CVX Toolbox and did a direct comparison with the Chebfun Remez algorithm (disclaimer: I am part of the Chebfun development team):

% Do the convex optimization bit.
m = 101; n = 11;            % 101 points, polynomial of degree 10
xi = linspace(-1, 1, m);    % equidistant points in [-1, 1]
ri = 1 ./ (1+(5*xi).^2);    % Runge function

tic                         % p is the polynomial of degree (n-1)
cvx_begin                   % minimize the distance in all points
variable p(n);
minimize( max(abs(polyval(p, xi) - ri)) );
cvx_end
toc_or = toc                % 0.17 sec for Matlab, CVX and SeDuMi

% Extract a Chebfun from the result
x = chebfun( [-1,1] );
A = [ chebfun(1) , x ];
for k=3:n, A(:,k) = A(:,k-1).*x; end
or = A * flipud(p)

% Make a chebfun of Runge's function
f = chebfun( @(x) 1 ./ ( 1 + 25*x.^2 ) )

tic, cr = remez( f , 10 ); toc_cr = toc

% Get the maximum error in each case
fprintf( 'maximum error of convex optimization: %e (%f s)\n' , norm( f - or , inf ) , toc_or );
fprintf( 'maximum error of chebfun remez: %e (%f s)\n' , norm( f - cr , inf ) , toc_cr );

% Plot the two error curves
plot( [ f - cr , f - or ] );
legend( 'chebfun remez' , 'convex optimization' );


After a lot of output, I get, on my laptop with Matlab 2012a, CVX version 1.22 and Chebfun's latest SVN Snapshot:

maximum error of convex optimization: 6.665479e-02 (0.138933 s)
maximum error of chebfun remez: 6.592293e-02 (0.309443 s)


Note that the Chebfun f used to measure the error is accurate to 15 digits. So Chebfun's Remez takes twice as long, but gets a smaller global error. It should be pointed out that CVX uses compiled code for the optimization whereas Chebfun is 100% native Matlab. The minimum error of 0.00654 is the minimum error 'on the grid', off that grid, the error can be up to 0.00659. Increasing the grid size to m = 1001 I get

maximum error of convex optimization: 6.594361e-02 (0.272887 s)
maximum error of chebfun remez: 6.592293e-02 (0.319717 s)


i.e. almost the same speed, but the discrete optimization is still worse as of the fourth decimal digit. Finally, ncreasing the grid size further to m = 10001 I get

maximum error of convex optimization: 6.592300e-02 (5.177657 s)
maximum error of chebfun remez: 6.592293e-02 (0.312316 s)


i.e. the discrete optimization is now more than ten times slower and is still worse as of the sixth digit.

The bottom line is that Remez will get you the globally optimal result. While the discrete analog may be fast on small grids, it will not give a correct result.

• And N. Trefethen was stressing the same and gave a similar example in the article I was citing. The question was not about the best approximation, but : What is the advantage of the Remez algorithm (nowadays) if you can get the same result with a reasonable convex solver? Mar 6 '12 at 7:55
• @HansWerner, I'm sorry, I misread your question. Your convex solver is not giving you the same result, at least not to all digits. If I understand your convex code correctly, you're minimizing the maximum error over a discrete set of points. This is a good approximation of -- but not the same as -- minimizing the global maximum error. Mar 6 '12 at 10:14
• The convex solver in this case gave a better result. Think of it, the Remez algorithm is an iterative procedure, quite similar to an optimization routine, and will not return an exact result, too. In the concrete case above, the solution from the optimization was better, i.e. had an overall smaller maximum norm, than the result from the best Remez implementation I know. The question is still open. Mar 6 '12 at 16:59
• @HansWerner, how did you measure the maximum error of the solution obtained with the convex solver? The Remez algorithm in chebfun should iterate until the minimum is achieved to machine precision (in a sense). Mar 6 '12 at 17:37
• Not necessarily; there are options like 'tol' (which is relative tolerance) or 'maxiter' for chebfun/remez, but there are similar options for all kinds of optimization solvers. In a way one could say, Remez is an optimization routine specialized for a certain task. And the question is: Have the general purpose solvers caught up and there is no need anymore for such a specialized solver? Of course, I may be wrong. Mar 6 '12 at 20:16