[R] trigonometric regression
William Dunlap
wdunlap at tibco.com
Thu Jun 17 18:12:07 CEST 2010
> -----Original Message-----
> From: r-help-bounces at r-project.org
> [mailto:r-help-bounces at r-project.org] On Behalf Of Duncan Murdoch
> Sent: Thursday, June 17, 2010 3:19 AM
> To: William Simpson
> Cc: r-help at r-project.org
> Subject: Re: [R] trigonometric regression
>
> William Simpson wrote:
> > Suppose I do a trigonometric regression
> > fit<-lm(y~ cf + sf)
> > where cf and sf are the cos and sine components.
> >
> > b<-coef(fit)
> > I have the fitted sine component b[2] and the cos component b[3].
> > Doing summary(fit) gives me the p-values and SEs for b[2] and b[3].
> >
> > But I want the amplitude of the fitted waveform
> > amp<-sqrt(b[2]^2+b[3]^2)
> >
> > Can someone please tell me how to get the p-value for amp?
> >
>
> "the p-value for amp" is ambiguous; p-values refer to tests, not
> functions. But let's assume you want to test whether amp = 0. Then
> this is achieved by an F test comparing the model with cf and
> sf versus
> one without it. You'll see it in summary(fit) at the bottom of the
> display. If you want to include other covariates in the
> model, you can
> use anova, e.g.
>
> anova(lm(y ~ other), lm(y ~ cf + sf + other))
You can also define a function that keeps the cos
and sin terms together so anova(fit) shows
one entry for the (cos,sin) pair. E.g., define
the following function
cs <- function(x, freq)cbind(cos=cos(x*freq), sin=sin(x*freq))
and use it as in
> time<- sort(runif(30,0,20))
>
y<-sin(1.7*time+.1)+2*sin(1.1*time+.8)+rnorm(length(time),mean=2,sd=0.3)
> fit <- lm(y~cs(time,1.7)+cs(time,1.1)+cs(time,1.3))
> anova(fit)
Analysis of Variance Table
Response: y
Df Sum Sq Mean Sq F value Pr(>F)
cs(time, 1.7) 2 12.000 5.9998 50.7267 3.692e-09 ***
cs(time, 1.1) 2 57.355 28.6776 242.4609 3.493e-16 ***
cs(time, 1.3) 2 0.780 0.3902 3.2988 0.05501 .
Residuals 23 2.720 0.1183
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com
>
> Duncan Murdoch
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
More information about the R-help
mailing list