/**************     BOX-COX     ************************
Author: Marc Nerlove
			AREC/2200 Symons Hall
			University of Maryland
			College Park, Maryland 20742-5535
			301-405-1388
			FAX: 301-314-9032
			mnerlove@arec.umd.edu
			http://www.arec.umd.edu/mnerlove/mnerlove
Submitted: Nov 1998
Code provided gratis without performance guarantees for public
non-commercial use.
***********************************************************/
/*GENERAL PROCEDURE FOR THE BOX-COX TRANSFORM.*/
/*There is a serious problem with the Box-Cox transform near
zero on a finite precision computer.
The problem is that in the expression (x^a -  1)/a, x>0,
a<1, x^a reaches its limit of precision before the limit of
precision of division
by a. This produces what I have called a dorsal finned likelihood
function in
the case of the wine data where the parameter crosses zero. It
causes failure to
converge in the GAUSS library CML and problems with grid
search as well. A
naive replacement of the transform by LN( ) at zero simply will
not work
satisfactorily.*/
/*zz is any N by 1 vector. lambda is the Box-Cox parameter.*/
proc  boxcox( zz, lambda );
    local zbc;
    if abs(lambda)  .LE 1e-8;
zbc = ln(abs(zz));
    else;
zbc = ((zz^lambda)-1) ./ lambda;
    endif;
    retp(zbc);
endp;
/*********************End Box-Cox Procedure**********************/



/***************** supplementary code *************************
author: anonymous
date: 2 Nov 98
Code provided gratis without performance guarantees for public
non-commercial use.

The naive calculation of the Box-Cox function, as Professor Nerlove has
pointed out in a separate post, fails for lambda < 1e-8.  This might not
seem to be much of a problem, but it can play havoc with the calculation
of numerical derivatives when estimating models that include the Box-Cox
function using Maxlik or CML.   Professor Nerlove's solution, truncating
the
calculation at 1e-8, is problematic as well for calculating numerical
derivatives because it creates a discontinuity at 1e-8 and -1e-8.

I recommend the procedure below.  It is accurate to 14 places or so for
any value of lambda.
**************************************************************/

proc boxcox(x,lambda);
    for i(1,rows(x),1);
       for j(1,cols(x),1);
     x[i,j] = boxcox0(x[i,j],lambda);
      endfor;
    endfor;
endp;



 proc boxcox0( x, lambda );
   local a, b, result;

   if x == 1;
       retp(0);
   endif;

   if lambda > 1 or lambda < -1;
       retp( (exp(lambda*ln(x)) - 1) / lambda );
   elseif lambda == 1;
       retp( x - 1 );
   elseif lambda == 0;
       retp( ln(x) );
   else;

      a = ln(x);
      b = lambda * a;
      result = a;
      for i(2,100,1);
  a = a * b / i;
  if a == 0;
      break;
  endif;
  result = result + a;
      endfor;
      retp(result);

   endif;
endp;