How to fit a curve with an equation

Hi Community,

I need to fit a curve with this equation : (-k*t+To)+t^-0,5

"k" and "To" are two constant and "t" is my time.

I did this fit on Excel (Photo_1) thanks to the solver with the below conditions :

  • Sum of square (SSR) : I want the minimum
  • I want to find the better "k" and "To"

However I don't know how I can do this on R.

My aim : I want to fit my curve on R with this equation and find k and To.

For example, I did that just below but I have three problems (and maybe it's not the good solution):

  • It give me a linear regression (photo_2 )
  • My coefficient aren't free (I need to find a solution where my coefficients can change in function of the better result found by R).
  • My last line code doesn't works

Photo_2

t <- PLR2decroissance$pupil_timestamp_secondes      # these are my time points
k <- 2.5         # time constant
To <- 0.48        # growth constant
y <- (-k*t+To)+t^-0.5 # my equation

# visualise
par(mfrow = c(1, 2))
plot(t, y)    
lines(PLR2decroissance$pupil_timestamp_secondes , PLR2decroissance$diameter_3d_corrige ,  col="red")

Thank you very much for your help.

Yours sincerely,

GIOVANNANGELI Cyril

It would help if you could post your data in a readable form, perhaps using dput(). Past that, you probably want to use lm() to estimate k and To.

Hi Startz,

Thank you for your answer. When you sayed " It would help if you could post your data in a readable form, perhaps using dput() " : you want I share my datas on the chat (if it's that, I can try).

However I'm novice on R and I don't understand how I can estimate k and To. I'm lost ..

Thank you for your help.

Yours sincerely,

Try copying an pasting your data in reply. Make sure that it's something that the reader can copy and paste; in other words, not picture.

I don't know how I can copying and pasting my datas (I'm sorry but I'm verry new in this field). So I exported my dataframe in a csv and I copy and past 64 row of my 2 columns "pupil_timestamp_secondes" and "diameter_3d_corrige". Is it works like that ?

pupil_timestamp_secondes
870.485044999979
870.491905000003
870.499957999971
870.507964999997
870.515828999982
870.524032999994
870.531897999987
870.539931999985
870.547885000007
870.559999999998
870.567870999977
870.575948999962
870.583884999971
870.591989999986
870.59990500001
870.608161000011
870.615955999994
870.62399599998
870.631705999956
870.640095999988
870.647958999965
870.659928000008
870.667883999995
870.67597099999
870.683879999968
870.692008999991
870.700002999976
870.70792700001
870.715874999994
870.724040000001
870.731809999968
870.739983999985
870.747930000012
870.759995999979
870.767835999955
870.775861000002
870.783901999996
870.791981999995
870.799994000001
870.808032999979
870.816004999971
870.824086999986
870.831732999999
870.840151999961
870.847921999986
870.860103999963
870.86785499996
870.875912999967
870.883927999996
870.892128999985
870.899872999988
870.908027999976
870.915923999972
870.924155999965
870.932074000011
870.939974999987
870.947936999961
870.960197999957
870.968147999956
870.975954999973
870.983956000011
870.991959999956
870.999833000009
871.007992999977

diameter_3d_corrige
1.48164707707831
1.49456095678936
1.48414707057656
1.48726408870283
1.48205378247385
1.47781943384314
1.47231104221379
1.46745207409162
1.45597322337438
1.43931190797885
1.42335791963431
1.42252471948163
1.39743728364405
1.39648141401239
1.38202497575066
1.36839658349799
1.3533448389161
1.34086410176018
1.32458607292278
1.32177503309484
1.31662963954374
1.30383311815029
1.28886017800971
1.28617193356825
1.2679510253103
1.26079463527277
1.25409525242261
1.2433266606438
1.23210515620124
1.21748213073997
1.21730916830565
1.20922005151322
1.19459113838058
1.18862563974446
1.18891287365467
1.1827331916104
1.16622961119455
1.16675354050288
1.16004219806803
1.15531594345702
1.15201180284821
1.14496753650957
1.14352465118899
1.12535131292713
1.12543821423565
1.12961750763861
1.11209925382821
1.11601875932558
1.11436788570255
1.1097258568122
1.10941190041669
1.10137899634844
1.10550057001415
1.09148442833974
1.10311317989005
1.08919209487414
1.09063657884501
1.08955663942749
1.09088134950307
1.09373823583395
1.09270374068362
1.07977138687109
1.08426291442292
1.08552464145148

That helps.
Try
lm(diameter_3d_corrige-1/sqrt(pupil_timestamp_secondes) ~ pupil_timestamp_secondes)
to get values for k and To.

I made up raw data x, y

set.seed(20160227)
x<-seq(.01,1,.01)
y<-runif(1,5,100)exp(-runif(100,0.01,0.05)x)+rnorm(100,0,0.5)
plot(x,y, col="blue")
formula <- y~ (-a
x+b)+x^-0.5
a_start <- .48
b_start <- 2.57
m<-nls(y ~ I(a
exp(-b*x)+ x^-0.5),start=list(a=a_start,b=b_start))
lines(x,predict(m),col="red",lty=2,lwd=3)

Thank you very much for your answer. I was on the good road, because yesterday, I tried with the "nls" function but I don't understand how use this function, even if I use internet. So thank you for your help.

Now, I tried your solution since 30 min :

x = PLR2decroissance$pupil_timestamp_secondes
y = PLR2decroissance$diameter_3d_corrige

plot(x,y, col="blue")        # This part works for me
formula <- y~ (-ax+b)+x^-0.5
a_start <- .48
b_start <- 2.57
m<-nls(y ~ I(a exp(-b*x)+ x^-0.5),start=list(a=a_start,b=b_start))       # It appears an error for this row
lines(x,predict(m),col="red",lty=2,lwd=3)

Now I have this message : " Error: unexpected symbol in "m<-nls(y ~ I(a exp" "

In addition, sorry, but I don't understand why we use "a_start" and "b_start" (I saw that at this link nls function - RDocumentation but I don't understand the uses) and this part with the "I" function " I(a exp "

I'm thinking, after, when it works, I will have in my dataframe "m" the coefficients a and b ?

Let me point out that this is actually a linear regression. No need to use nls() to get the coefficient values.

Thank you startz, I tried this below and it works

However, if I want to use the exposed formula in the beginning, how I can do it ? Because I showed to my advisor some weeks ago a linear regression on this curve and he sayed in this case it's better to use a formula with a linear regression + an exponential because we have a curve at the ending of the graph.

You said you wanted

That's what you now have. You've estimated k and To exactly for that equation. (Well, actually the coefficient from the linear regression is -k, (my fault)).

Maybe it's my fault and I don't understand .. but " y <- (-k*t+To)+t^-0.5 " is not a linear equation.

Because when I do that in excel, I have a curve (Photo_4). The fiting in Excel is not very good but I would like do that with R

Photo_4

Currently, with your (Photo_6) formula and an other formula (Photo_5), I have two differents linear regression, so I have a straight line and not a curve :

Photo_5

Photo_6

Sorry for my ununderstanding and if my questions looks like lame

Your questions are perfectly good ones.

A regression is linear if it is linear in the parameters to be estimated. It doesn't matter if some of the things on the right happen to be nonlinear functions so long as the parameters to be estimated are outside those functions. That's the case in your model.

To get the parameters it turned out to be easiest to take the square root term to the left. That means that the predicted values you used in abline aren't exactly what you want.

By the way, are you sure you don't want a coefficient in front of t^-0.5 in your equation?

Thank you very much for these explanations. So if I understand, even if my equation look like non-linear, if my parameters are linear .. I will have a linear function. So now, I have an other question : why in excel I have a curve for this formula ?

Maybe a coefficient in front of t^-0.5 in my equation will be better ? What is your proposition ?

And about the abline function, I'm understanding, this function graph a straight line. So maybe it's the reason why it draw a straight line .. but with the "lines" function it doesn't works also, it draw a straight line also. But if I think about your explanations, it's normal because my parameters are linear. So I come back to my question on the top of this post, why in excel I have a curve for this formula ? :sweat_smile:

Thank your for you patience in you explanations.

I'm not sure what you did in Excel. I note that some of the pictures posted are for different data. The vertical scales are differernt.

You can try something like

lm(diameter_3d_corrige ~ I(1/sqrt(pupil_timestamp_secondes)) + pupil_timestamp_secondes)

Thank you, your code works but it give a straight line.

I would like realise a curve like this video : Nonlinear curve fitting in R using mosaic and nls - YouTube. It explain, we can fit a straight line but we can also find the curve who fitting with datas. I did what it say in this video, but it doesn't works an other once, I have a straight line.

The code is really similar to @fcas80 and I share to you what I did :

library(mosaic)

#create data frame
decroissance <- data.frame(x=PLR2decroissance$pupil_timestamp_secondes,
                 y=PLR2decroissance$diameter_3d_corrige)
attach(decroissance)
plot(PLR2decroissance$pupil_timestamp_secondes, PLR2decroissance$diameter_3d_corrige, pch=15)

x = PLR2decroissance$pupil_timestamp_secondes
y = PLR2decroissance$diameter_3d_corrige


k_start <- .48
To_start <- 2.57
fitcurve = nls(y ~ (-k*x+To)+x^0.5, start=list(k=k_start,To=To_start))
coef(fitcurve)    
confint(fitcurve, level = 0.9)
plot(x, y, pch=15)
lines(x, predict(fitcurve), col="blue")

It's not a straight line. It has the term x^0.5 in it. However, that term may be so small that it still looks like a straight line to the eye.

Note that you are missing a minus sign on the x^0.5

Also note that your scale in Photo_4 is very different from the data you provided.

I suspect that Excel was solving a different problem.

1 Like

Thank you very much for your answer. With a help, I found the code, so I share the code who can give me a curve like the curve on Excel, maybe you can understand me better :

# y ~ (-k*x+To)+x^-0.5

mod <- glm(y ~ pupil_timestamp_secondes + sqrt(pupil_timestamp_secondes), df, family="gaussian")
summary(mod)
plot(y~pupil_timestamp_secondes, df, type="l")
lines(df$pupil_timestamp_secondes, predict(mod), col = 'blue')

And you can see the result on the picture

Now, I would like do the same for this curve just below with this formula, but it's very more difficult :sweat_smile:

Formula : A1*(1-exp(-t/To1))+A2*(1-exp(-t/To2))+D0

Picture :

Someone can help me ? :sweat_smile:

Yours sincerely,

GIOVANNANGELI Cyril

I tried that just below but R send me this message "Error in nls ... singular gradient". Maybe someone have a correction ? (I share my code below)

If you want I can share my datas in csv with this wetransfert link : WeTransfer - Send Large Files & Share Photos Online - Up to 2GB Free

x = df5$pupil_timestamp_secondes
y = df5$y

A1 <- 3; lrc1 <- 1; A2 <- 400; lrc2 <- 4000
mod <- nls(y ~ SSbiexp(x, A1, lrc1, A2, lrc2), data = df5 , trace = TRUE)
summary(mod)
plot(y~x, df5, type="l")
lines(df5$x, predict(mod), col = 'blue')

Yours sincerely,

GIOVANNANGELI Cyril

Nobody have a solution ? :sweat_smile: