[R] Hodrick-Prescott Filter: Special Case of smooth.spline()

[Martin Maechler]

>>>>> "Roger" == Roger Koenker <roger at ysidro.econ.uiuc.edu> writes:

    Roger> A general answer would be to forget the
    Roger> Hodrick-Prescott numerology for lambda=1600, and use
    Roger> the spar parameter produced by smooth.spline's GCV
    Roger> functionality.  Nevertheless, it would be nice for
    Roger> reproducibility (comparison) reasons to be able to
    Roger> get the actual lambda corresponding to spar.  My
    Roger> quick look at the code didn't yield an easy way to
    Roger> extract the ratio: r = tr(X' W^2 X) / tr(Sigma), and
    Roger> thereby get back to lambda.  Maybe Brian could
    Roger> suggest something for this.

and Brian privately mentioned that I had most recently worked on
this.  Exactly because I found I did want to get at the lambda
used, in R 1.4 (to be released in a few days) `lambda' at least
is returned as well as spar and the details section contains
even more details about the relation of `spar' and `lambda'
((and `spar' is no longer restricted strictly to [0,1]; "0" particularly
  came from an S ``pseudo-compatibility'' which does not make
  much sense for our spar which is on a LOG scale)).

    Roger> I might also note that the Hodrick Prescott setup
    Roger> would require that you set all.knots=TRUE, since the
    Roger> default in smooth.spline is to chose fewer knots.
    Roger> The HP penalty assumes equally spaced x's so the
    Roger> penalty is just the sum of the second differences of
    Roger> ghat.

    Roger> For those of you (if any) who are wondering what the
    Roger> Hodrick-Prescott filter is: in effect it is the
    Roger> special case of the (Reinsch) cubic smoothing spline
    Roger> for an equally spaced time series.  HP claim that for
    Roger> (typical quarterly economic time-series) lambda can
    Roger> be chosen to be 1600, and (amazingly) this has become
    Roger> a default smoother in some circles of
    Roger> macroeconometrics.

interesting... this means that one uses the same (kernel-) equivalent
smoothing window for all quartely time series?  Wouldn't this
mean that (a combination of) the smoothness of the true
underlying function m(x) and the (local) variance sigma(x) was
assumed to be (almost) a universal constant?
interesting, indeed...

    Roger> url: http://www.econ.uiuc.edu Roger Koenker email
    Roger> roger at ysidro.econ.uiuc.edu Department of Economics
    Roger> vox: 217-333-4558 University of Illinois fax:
    Roger> 217-244-6678 Champaign, IL 61820

    Roger> On Fri, 14 Dec 2001, Nick Davis wrote:

    >> I've had a play with this and, due to my own
    >> short-comings, remain none the wiser.
    >> 
    >> In particular, I'm not sure what value of 'spar' is
    >> consistent with the magic lambda=1/1600 for quarterly
    >> data.
    >> 
    >> I initially interpreted spar as lambda and tried setting
    >> spar=1/1600.  This results in almost no smoothing while
    >> spar=1600 causes an error.  The smooth.spline function
    >> seems to want something in the region of 0 to 1.  The
    >> closer spar is to 1, the greater the smoothing.  A closer
    >> look at the R documentation revealed that:
    >> 
    >> The computational lambda used (as a function of `spar')
    >> is lambda = r * 256^(3*spar - 1) where r = tr(X' W^2 X) /
    >> tr(Sigma), Sigma is the matrix given by Sigma[i,j] =
    >> Integral B''[i](t) B''[j](t) dt, X is given by X[i,j] =
    >> B[j](x[i]), W^2 is the diagonal matrix of scaled weights,
    >> `W = diag(w)/n' (i.e., the identity for default weights),
    >> and B[k](.) is the k-th B-spline.
    >> 
    >> I admit to being a bit confused by the matrix algebra.
    >> It appears to come down to knowing 'r' so that 'spar' can
    >> be derived by imposing a constraint on lambda.

yes;  as `r' only depends on the design (i.e. x[] and weights[]),
in R 1.4 you can use one value of spar, and from the result of
smooth.spline compute r from spar and lambda.
As Brian mentioned (in a private mail), we could allow `lambda'
as an alternative argument to `spar' directly.  
Not for 1.4 though since that's in deep feature freeze.


    >> If anyone can shed some light on this it would be much
    >> appreciated.  A general answer would be nice as I don't
    >> always work with quarterly data.
    >> 
    >> Thanks & Regards,
    >> 
    >> Nick Davis


  < .............. >


Martin Maechler <maechler at stat.math.ethz.ch>	http://stat.ethz.ch/~maechler/
Seminar fuer Statistik, ETH-Zentrum  LEO C16	Leonhardstr. 27
ETH (Federal Inst. Technology)	8092 Zurich	SWITZERLAND
phone: x-41-1-632-3408		fax: ...-1228			<><
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._