[Rd] poly() can exceed degree k - 1 for k distinct points (PR#11251)
Steven McKinney
smckinney at bccrc.ca
Fri Apr 25 05:26:55 CEST 2008
I haven't parsed the source to fully understand
the 'normalization constants' returned in the
poly output component
attr(,"coefs")$norm2
but notice that the first 6 are non-zero and
the last 8 are smaller than machine precision.
Some kind of useful warning to the user could
ensue upon evaluation of such normalization
constants?
> x = rep(1:5,3)
> x
[1] 1 2 3 4 5 1 2 3 4 5 1 2 3 4 5
> poly(x, 12)
1 2 3 4 5 6 7 8 9 10
[1,] -3.651484e-01 0.3086067 -1.825742e-01 0.06900656 0.50410540 -0.04690012 0.16698778 0.03700268 -0.21853454 -0.50410540
[2,] -1.825742e-01 -0.1543033 3.651484e-01 -0.27602622 0.11156408 -0.48312338 -0.10874588 -0.02250863 0.27412822 -0.11156408
[3,] -6.437201e-18 -0.3086067 -1.301129e-16 0.41403934 -0.17911692 -0.14493516 0.51592172 0.11434917 0.03151684 0.17911692
[4,] 1.825742e-01 -0.1543033 -3.651484e-01 -0.27602622 0.07326369 0.11056157 -0.07141299 0.52748012 0.18001893 -0.07326369
[5,] 3.651484e-01 0.3086067 1.825742e-01 0.06900656 -0.17121360 -0.25395954 -0.14949152 0.20036589 -0.42074883 0.17121360
[6,] -3.651484e-01 0.3086067 -1.825742e-01 0.06900656 -0.75205270 0.02345006 -0.08349389 -0.01850134 0.10926727 -0.24794730
[7,] -1.825742e-01 -0.1543033 3.651484e-01 -0.27602622 -0.05578204 0.74156169 0.05437294 0.01125431 -0.13706411 0.05578204
[8,] -6.437201e-18 -0.3086067 -8.788409e-17 0.41403934 0.08955846 0.07246758 -0.75796086 -0.05717459 -0.01575842 -0.08955846
[9,] 1.825742e-01 -0.1543033 -3.651484e-01 -0.27602622 -0.03663185 -0.05528079 0.03570650 -0.76374006 -0.09000946 0.03663185
[10,] 3.651484e-01 0.3086067 1.825742e-01 0.06900656 0.08560680 0.12697977 0.07474576 -0.10018295 0.71037442 -0.08560680
[11,] -3.651484e-01 0.3086067 -1.825742e-01 0.06900656 0.24794730 0.02345006 -0.08349389 -0.01850134 0.10926727 0.75205270
[12,] -1.825742e-01 -0.1543033 3.651484e-01 -0.27602622 -0.05578204 -0.25843831 0.05437294 0.01125431 -0.13706411 0.05578204
[13,] -6.437201e-18 -0.3086067 -8.788409e-17 0.41403934 0.08955846 0.07246758 0.24203914 -0.05717459 -0.01575842 -0.08955846
[14,] 1.825742e-01 -0.1543033 -3.651484e-01 -0.27602622 -0.03663185 -0.05528079 0.03570650 0.23625994 -0.09000946 0.03663185
[15,] 3.651484e-01 0.3086067 1.825742e-01 0.06900656 0.08560680 0.12697977 0.07474576 -0.10018295 -0.28962558 -0.08560680
11 12
[1,] 0.04690012 0.16698778
[2,] 0.48312338 -0.10874588
[3,] 0.14493516 0.51592172
[4,] -0.11056157 -0.07141299
[5,] 0.25395954 -0.14949152
[6,] -0.02345006 -0.08349389
[7,] 0.25843831 0.05437294
[8,] -0.07246758 0.24203914
[9,] 0.05528079 0.03570650
[10,] -0.12697977 0.07474576
[11,] -0.02345006 -0.08349389
[12,] -0.74156169 0.05437294
[13,] -0.07246758 -0.75796086
[14,] 0.05528079 0.03570650
[15,] -0.12697977 0.07474576
attr(,"degree")
[1] 1 2 3 4 5 6 7 8 9 10 11 12
attr(,"coefs")
attr(,"coefs")$alpha
[1] 3.000000 3.000000 3.000000 3.000000 3.000000 1.314957 2.355111 2.973300 4.032925 4.323708 1.314957 2.355111
attr(,"coefs")$norm2
[1] 1.000000e+00 1.500000e+01 3.000000e+01 4.200000e+01 4.320000e+01 2.468571e+01 1.128824e-28 1.344243e-28 7.418009e-29 9.824087e-27
[11] 1.017082e-32 2.357175e-51 9.543887e-57 6.351670e-57
attr(,"class")
[1] "poly" "matrix"
Steven McKinney
-----Original Message-----
From: r-devel-bounces at r-project.org on behalf of John Maindonald
Sent: Thu 4/24/2008 5:24 PM
To: r-devel at r-project.org
Subject: Re: [Rd] R-devel Digest, Vol 62, Issue 24
Actually, this may be a useful feature! It allows calculation of a
basis for the orthogonal complement of the space spanned by
model.matrix(lm(y ~ poly(x,12)). However, the default ought surely to
be to disallow df > k-1 in poly(x,df), where k = length(unique(x)).
John Maindonald email: john.maindonald at anu.edu.au
phone : +61 2 (6125)3473 fax : +61 2(6125)5549
Centre for Mathematics & Its Applications, Room 1194,
John Dedman Mathematical Sciences Building (Building 27)
Australian National University, Canberra ACT 0200.
On 24 Apr 2008, at 8:00 PM, r-devel-request at r-project.org wrote:
> From: russell-lenth at uiowa.edu
> Date: 24 April 2008 3:05:28 AM
> To: r-devel at stat.math.ethz.ch
> Cc: R-bugs at biostat.ku.dk
> Subject: [Rd] poly() can exceed degree k - 1 for k distinct points
> (PR#11251)
>
>
> The poly() function can create more variables than can be fitted when
> there are replicated values. In the example below, 'x' has only 5
> distinct values, but I can apparently fit a 12th-degree polynomial
> with
> no error messages or even nonzero coefficients:
>
> R> x = rep(1:5,3)
> R> y = rnorm(15)
> R> lm(y ~ poly(x, 12))
>
> Call:
> lm(formula = y ~ poly(x, 12))
>
> Coefficients:
> (Intercept) poly(x, 12)1 poly(x, 12)2 poly(x, 12)3
> -0.27442 0.35822 -0.26412 2.11780
> poly(x, 12)4 poly(x, 12)5 poly(x, 12)6 poly(x, 12)7
> 1.83117 -0.09260 -0.48572 1.94030
> poly(x, 12)8 poly(x, 12)9 poly(x, 12)10 poly(x, 12)11
> -0.88297 -1.04556 0.74289 -0.01422
> poly(x, 12)12
> -0.46548
>
<snip>
<snip>
> [I thought I submitted this via the website yesterday, but I can
> find no
> trace of it. I apologize if this is a duplicate, but I don't think
> it is.]
> --
> Russell V. Lenth, Professor
> Department of Statistics
> & Actuarial Science (319)335-0814 FAX (319)335-3017
> The University of Iowa russell-lenth at uiowa.edu
> Iowa City, IA 52242 USA http://www.stat.uiowa.edu/~rlenth/
[[alternative HTML version deleted]]
______________________________________________
R-devel at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
More information about the R-devel
mailing list