[R] fitting cosine curve

Charles C. Berry ccberry at ucsd.edu
Wed Jun 21 18:48:41 CEST 2017


On Wed, 21 Jun 2017, J C Nash wrote:

> Using a more stable nonlinear modeling tool will also help, but key is to get
> the periodicity right.
>

The model is linear up to omega after transformation as Don and I noted.

Taking a guess that 2*pi/240 = 0.0262 is about right for omega:

> rsq <- function(x) {t2<-t*x;summary(lm(y~cos(t2)+sin(t2)))$r.squared}
> vrsq <- Vectorize(rsq)
> optimise(rsq, c(0.8,1.2)*2*pi/240,maximum=TRUE)
$maximum
[1] 0.02794878

$objective
[1] 0.8127072

> curve(vrsq,0.025,0.03)
>

Isn't this stable enough?

And as you note

plot(lm(y~cos(t*0.0279)+sin(t*0.0279)))

reveals lack-of-fit.

Of course there are some other issues not addressed here such as 
possible autoregression.

HTH,

Chuck

> y=c(16.82, 16.72, 16.63, 16.47, 16.84, 16.25, 16.15, 16.83, 17.41, 17.67,
> 17.62, 17.81, 17.91, 17.85, 17.70, 17.67, 17.45, 17.58, 16.99, 17.10)
> t=c(7,  37,  58,  79,  96, 110, 114, 127, 146, 156, 161, 169, 176, 182,
> 190, 197, 209, 218, 232, 240)
>
> lidata <- data.frame(y=y, t=t)
>
> #I use the method to fit a curve, but it is different from the real curve,
> #which can be seen in the figure.
> linFit  <- lm(y ~ cos(t))
> library(nlsr)
> #fullFit <- nls(y ~ A*cos(omega*t+C) + B,
> #start=list(A=coef(linFit)[1],B=coef(linFit)[2],C=0,omega=.4))
> #omega cannot be set to 1, don't know why.
> fullFit <- nlxb(y ~ A*cos(omega*t+C) + B, data=lidata,
> start=list(A=coef(linFit)[1],B=coef(linFit)[2],C=0,omega=.04), trace=TRUE)
> co <- coef(fullFit)
> fit <- function(x, a, b, c, d) {a*cos(b*x+c)+d}
> plot(x=t, y=y)
> curve(fit(x, a=co['A'], b=co['omega'], c=co['C'],d=co['B']), add=TRUE
> ,lwd=2, col="steelblue")
> jstart <- list(A=20, B=100, C=0, omega=0.01)
> jfit <- nlxb(y ~ A*cos(omega*t+C) + B, data=lidata,
>        start=jstart, trace=TRUE)
> co <- coef(jfit)
> fit <- function(x, a, b, c, d) {a*cos(b*x+c)+d}
> plot(x=t, y=y)
> curve(fit(x, a=co['A'], b=co['omega'], c=co['C'],d=co['B']), add=TRUE
>      ,lwd=2, col="steelblue")
>
> JN
>
> On 2017-06-21 12:06 AM, lily li wrote:
>> I'm trying the different parameters, but don't know what the error is:
>> Error in nlsModel(formula, mf, start, wts) :
>>    singular gradient matrix at initial parameter estimates
>> 
>> Thanks for any suggestions.
>> 
>> On Tue, Jun 20, 2017 at 7:37 PM, Don Cohen <don-r-help at isis.cs3-inc.com>
>> wrote:
>> 
>>> 
>>> If you know the period and want to fit phase and amplitude, this is
>>> equivalent to fitting a * sin + b * cos
>>>
>>>   > >>> > I don't know how to set the approximate starting values.
>>> 
>>> I'm not sure what you meant by that, but I suspect it's related to
>>> phase and amplitude.
>>>
>>>   > >>> > Besides, does the method work for sine curve as well?
>>> 
>>> sin is the same as cos with a different phase
>>> Any combination of a and b above = c * sin (theta + d) for
>>> some value of c and d and = e * cos (theta + f) for some value
>>> of e and f.
>>> Also for any c,d and for any e,f there is an a,b.
>>> the c and e are what I'm calling amplitude, the d and f are what
>>> I'm calling phase.
>>> 
>>
>> 	[[alternative HTML version deleted]]
>> 
>> ______________________________________________
>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> 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.
>> 
>

Charles C. Berry                 Dept of Family Medicine & Public Health
cberry at ucsd edu               UC San Diego / La Jolla, CA 92093-0901
http://biostat.ucsd.edu/ccberry.htm



More information about the R-help mailing list