[R] Specify model with polynomial interaction terms up to degree n
Rui Barradas
ruipbarradas at sapo.pt
Mon Jul 2 19:14:30 CEST 2012
Hello,
Another way is to cbind the vectors 'a' and 'b', but this needs argument
'raw' set to TRUE.
poly(cbind(a, b), 6, raw=TRUE)
To the OP: is this time series related? With 6 being a lag or test
(e.g., Tsay, 1986) order? I'm asking this because package nlts has a
function for this test up to order 5 and it uses poly().
Hope this helps,
Rui Barradas
Em 02-07-2012 16:04, David Winsemius escreveu:
>
> On Jul 2, 2012, at 10:51 AM, David Winsemius wrote:
>
>>
>> On Jul 2, 2012, at 9:29 AM, YTP wrote:
>>
>>> I would like to specify a model with all polynomial interaction terms
>>> between
>>> two variables, say, up to degree 6. For example, terms like a^6 + (a^5 *
>>> b^1) + (a^4 * b^2) + ... and so on. The documentation states
>>>
>>> The ^ operator indicates crossing to the specified degree.
>>>
>>> so I would expect a model specified as y ~ (a+b)^6 to produce these
>>> terms.
>>> However doing this only returns four slope coefficients, for
>>> Intercept, a,
>>> b, and a:b. Does anyone know how to produce the desired result?
>>> Thanks in
>>> advance.
>>
>> You might try:
>>
>> poly(a,6)*poly(b,6)
>>
>> (untested ... and it looks somewhat dangerous to me.)
>
> Well, now it's tested and succeeds at least numerically. Also tested
>
> ( poly(a,6) +poly(b,6) )^2 with identical results.
>
> Whether this is wise practice remains in doubt:
>
> dfrm <- data.frame(out=rnorm(100), a=rnorm(100), b=rnorm(100) )
> anova(lm( out ~ (poly(a,6) +poly(b,6) )^2, data=dfrm) )
> #-----------------------
> Analysis of Variance Table
>
> Response: out
> Df Sum Sq Mean Sq F value Pr(>F)
> poly(a, 6) 6 12.409 2.06810 3.0754 0.01202 *
> poly(b, 6) 6 5.321 0.88675 1.3187 0.26596
> poly(a, 6):poly(b, 6) 36 41.091 1.14142 1.6974 0.04069 *
> Residuals 51 34.295 0.67246
> ---
> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
More information about the R-help
mailing list