[R] Hodrick-Prescott Filter: Special Case of smooth.spline()
Martin Maechler
maechler at stat.math.ethz.ch
Fri Dec 14 09:54:11 CET 2001
>>>>> "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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
More information about the R-help
mailing list