How to calculate the coefficient "a" in a non-linear regression "y=2/(1 exp(a*x))-1&q-CodePudding

As stated in the topic, given some value of x (FTSW), and its corresponding value of y (NLE), with the euqation y=2/(1 exp(a*x))-1, how could I calcualate the "a" and draw the graph?

structure(list(FTSW = c(0.900359288405312, 0.900334904035704, 
                        0.905841413662097, 0.667621859862976, 0.632473107345453, 0.462837958498518, 
                        0.400045032939416, 0.352560790392534, 0.377856799586057, 0.628289954726284, 
                        0.311242672407523, 0.170933345859364, 0.311242672407523, 0.315689846065931, 
                        0.14105342778594, 0.147112396871317, 0.0722352096264757, 0.0615716495856242, 
                        0.396997072279772, 0.0332444780173884, 0.0871102539376406, 0.320745031679664, 
                        0.232853736349989, 0.195548432729584), NLE = c(0.557811005824568, 
                                                                       0.947807189426057, 0.676214815603592, 0.444390402075227, 0.97602256699577, 
                                                                       0.728254923413566, 0.306279885240675, 0.970762905436273, 0.611065629690272, 
                                                                       0.253424268258028, 0.557708161582854, 0.108086143843963, 0.727886237647626, 
                                                                       0.53056603773585, 0.54236648424231, 0.338367448826758, 0.111826531578449, 
                                                                       0.231650365670985, 0.409424366533147, 0.047092400998931, 0.19369651429764, 
                                                                       0.342990129203114, 0.288499226405363, 0.295346695941639)), class = c("tbl_df", 
                                                                                                                                            "tbl", "data.frame"), row.names = c(NA, -24L))

I tried:

  df <- read_excel("##.xlsx")
  x <- FTSW
  y <- NLE

  ggplot(df, aes(x, y)) 
  geom_point(pch=19, cex=3, ,col="Blue3")  
  geom_smooth(method="nls", formula=y ~ 2/(1 exp(a*x))-1, method.args=list(start = list(a=1)), se=FALSE,lwd= 1.2, col= "black") 
  scale_y_continuous(limits = c(0,1.0))  
  scale_x_continuous(limits = c(0,1))  
  labs(title ="Blancs",y="Expansion folliaire verte relative",x="FTSW") 
  theme_update(plot.title=element_text(hjust=0.5))

With this, I can get the graph, but I don't know how to have the "a"value. Could anyone help me?

CodePudding user response：

You can estimate the model (say, mod) separately, using stat::nls, and obtain the coefficient via coef(mod)

mod = nls(NLE ~ 2/(1 exp(a*FTSW))-1,start = list(a=1),data = dat)
dat$pred = predict(mod,dat)
a = coef(mod)

ggplot(dat,aes(x=FTSW))   
  geom_point(aes(y=NLE), pch=19, cex=3, ,col="Blue3")  
  geom_line(aes(y=pred), color="black", lwd=1.2)   
  labs(caption=paste("'a' value from model=",round(a,4)))   
  scale_y_continuous(limits = c(0,1.0))  
  scale_x_continuous(limits = c(0,1))  
  labs(title ="Blancs",y="Expansion folliaire verte relative",x="FTSW") 
  theme_update(plot.title=element_text(hjust=0.5))

CodePudding user response：

Add a trace to nls which stores the result of running it in .nls in the global environment.

trace(nls, exit = quote(.nls <<- returnValue()))
# run ggplot
untrace(nls)

coef(.nls)
##         a 
## -2.974184