I have a continuous response variable (leaf length), which has been repeatedly measured in 45 experimental plots over several growing seasons. My goal is to properly estimate at which day of the year (DOY) certain events happened in each year (e.g., the passing of 20% or 80% of annual leaf growth). I fitted different smoothers for each plot using the gam function of mgcv. Here some examples for one year of 16 plots (x-axis: DOY; y-axis: leaf length in percent):
I then extracted the DOY and its confidence interval (CI) for the passing of each threshold by sampling the posterior according to this post: Can I use bootstrapping to estimate the uncertainty in a maximum value of a GAM?
So far so good. I have an estimate for the mean DOY with a lower and upper boundary of the CI for each plot and threshold. But now I would like to fit a mixed effect model to assess whether the timing of these events differs between treatments. I was looking for a way to take the uncertainty coming from each plot-level estimate into account.
My question: Is it correct to calculate the inverse of the squared standard error (derived from the CI) as a weight and use this with lme as following?
lme1 <- lme(DOY ~ year * treatment, random = ~1|plot, weights = varFixed(~weight), data = leaf_growth)
It plays a large role whether or not these weights are included.
OR is there an entirely different way to get good estimates for the passing of such thresholds?
