Re: [Rd] Rd] Numerics behind splineDesign

From: Nathaniel Smith <njs_at_pobox.com>
Date: Thu, 02 Aug 2012 15:11:59 +0100

On Thu, Aug 2, 2012 at 2:09 PM, Terry Therneau <therneau_at_mayo.edu> wrote:
>
>
> On 08/02/2012 05:00 AM, r-devel-request@r-project.org wrote:
>>
>> Now I just have to grovel over the R code in ns() and bs() to figure
>> out how exactly they pick knots and handle boundary conditions, plus
>> there is some code that I don't understand in ns() that uses qr() to
>> postprocess the output from spline.des. I assume this is involved
>> somehow in imposing the boundary conditions...
>>
>> Thanks again everyone for your help,
>> -- Nathaniel
>
> The ns and bs function post-process the spline bases to get an orthagonal
> basis matrix, this is the use of qr. I think this causes much more grief
> than it is worth, for the sake of a small increase in numeric stability.
> For instance when you plot the spline bases, they don't look anything like
> the basis functions one would expect. (Perhaps my background in numerical
> analysis was a hindrance here, since I know something about splines and thus
> have an expectation).

You know, the white book does claim that S orthogonalizes the spline basis functions, but R doesn't seem to actually do that:

> crossprod(bs(1:10, 4, intercept=TRUE))

           1 2 3 4

1 1.84104162 0.6154211 0.2570069 0.06430817
2 0.61542109 0.7710207 0.5787736 0.25700689
3 0.25700689 0.5787736 0.7710207 0.61542109
4 0.06430817 0.2570069 0.6154211 1.84104162

> crossprod(ns(1:10, 4, intercept=TRUE))
1 2 3 4 1 1.0659929 0.5314089 0.37125505 -0.17838443 2 0.5314089 1.0681360 0.23122018 0.14422236
3 0.3712550 0.2312202 1.30226657 -0.03809082 4 -0.1783844 0.1442224 -0.03809082 1.10144173

and they look quite reasonable when plotted to me.

ns is doing something more complicated; I think it's computing the second derivatives of each basis function at the position of the boundary knots, and then using that to somehow transform the original basis. Since the key feature of ns is that the returned basis has a zero second derivative at the boundary knots, I'm sure this makes sense if one stares at it for long enough.

-n



R-devel_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel Received on Thu 02 Aug 2012 - 14:15:06 GMT

This quarter's messages: by month, or sorted: [ by date ] [ by thread ] [ by subject ] [ by author ]

All messages

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Thu 02 Aug 2012 - 17:30:35 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-devel. Please read the posting guide before posting to the list.

list of date sections of archive