How to fit sigmoidal curves and extract x at a specified value of y using Bayesian statistics?

47 Views Asked by At

I have a dataset of multiple trees with measurement values between 0 and 1 taken over a multiple days. It is known that when plotted as day (x) and value (y), each tree's measurements form a sigmoidal curve (see graph). I have the following code to fit the curve using non-linear least squares, and extract the x value (day) when y passes through 0.5 (and the gradient of the curve at this point) as shown on the graph. We will then use these x values as a response variable in a spatial model, so we want to produce them using Bayesian methods instead, to account for the uncertainty in their estimation. enter image description here

I am a complete beginner at Bayesian so would like some help on how to modify this code so the model is produced with brms? I can't find any help pages showing how to do this that I understand.

Code:

results<-data.frame()
#iterate through each TreeID
trees<-unique(all$TreeID)
for (tree in trees) {
  #subset data for current TreeID
  subset_data<-filter(all, TreeID==tree)
  #fit sigmoidal logistic curve for current TreeID
  model<-nls(Value~SSlogis(Day, Asym, xmid, scal), data=subset_data)
  #extract parameters for current TreeID
  params<-coef(model)
  #calculate x when y=0.5 for current TreeID
  bbday<-params["xmid"]+params["scal"]*log(0.5/(params["Asym"]-0.5))
  #calculate gradient (derivative) when y=0.5 for current TreeID
  bbgrad<-params["Asym"]*exp((params["xmid"]-bbday)/params["scal"])/(params["scal"]*(exp((params["xmid"]-bbday)/params["scal"])+1)^2) 
  #store results in dataframe
  results<-rbind(results, data.frame(TreeID=tree, BudBurstDay=bbday, BudBurstGrad=bbgrad))
}

Here's data for the tree in the graph:

TreeID<-c(14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511, 14511)
Day<-c(84, 87, 89, 93, 96, 98, 103, 105, 108, 111, 114, 117, 120, 123, 127, 129, 132, 134, 138, 141, 144, 147, 150, 155, 158)
Value<-c(0.04166667, 0.06250000, 0.16666667, 0.26666667, 0.26666667, 0.30000000, 0.36666667, 0.42857143, 0.16666667, 0.43333333, 0.50000000, 0.41666667, 0.61111111, 0.77083333, 0.89583333, 1.00000000, 1.00000000, 1.00000000, NA, 1.00000000, 1.00000000, 1.00000000, 1.00000000, NA, 1.00000000)
0

There are 0 best solutions below