Thank you. This topic is more complicated than anticipated. Best regards, Luigi
On Tue, Apr 15, 2025 at 11:09 PM Andrew Robinson <a...@unimelb.edu.au> wrote: > > A statistical (off-topic!) point to consider: when the GLM was fitted, you > conditioned on x and let y be the random variable. Therefore the model > supports predictions of y conditional on x. You’re now seeking to make > predictions of x conditional on y. That’s not the same thing, even in OLS. > > It might not matter for your application but it’s probably worth thinking > about. Simulation experiments might shed some light on that. > > Cheers, Andrew > > -- > Andrew Robinson > Chief Executive Officer, CEBRA and Professor of Biosecurity, > School/s of BioSciences and Mathematics & Statistics > University of Melbourne, VIC 3010 Australia > Tel: (+61) 0403 138 955 > Email: a...@unimelb.edu.au > Website: https://researchers.ms.unimelb.edu.au/~apro@unimelb/ > > I acknowledge the Traditional Owners of the land I inhabit, and pay my > respects to their Elders. > > On 16 Apr 2025 at 1:01 AM +1000, Gregg Powell via R-help > <r-help@r-project.org>, wrote: > > > Take a look at this Luigi... > > > > # The model is: logit(p) = β₀ + β₁*Cycles > # Where p is the probability (or normalized value in your case) > > # The inverse function would be: Cycles = (logit⁻¹(p) - β₀)/β₁ > # Where logit⁻¹ is the inverse logit function (also called the expit > >function) > > # Extract coefficients from your model > intercept <- coef(b_model)[1] > slope <- coef(b_model)[2] > > # Define the inverse function > inverse_predict <- function(p) { > # p is the probability or normalized value you want to find the >cycles for > # Inverse logit: log(p/(1-p)) which is the logit function > logit_p <- log(p/(1-p)) > > > > > > > # Solve for Cycles: (logit(p) - intercept)/slope > cycles <- (logit_p - intercept)/slope > > > > > > > return(cycles) > } > > # Example: What cycle would give a normalized value of 0.5? > inverse_predict(0.5) > > > > > This function takes a probability (normalized value between 0 and 1) and > returns the cycle value that would produce this probability according to your > model. > Also: > This works because GLM with binomial family uses the logit link function by > default > The inverse function will return values outside your original data range if > needed > This won't work for p=0 or p=1 exactly (you'd get -Inf or Inf), so you might > want to add checks > > All the best, > Gregg > > > > > > > > > > On Tuesday, April 15th, 2025 at 5:57 AM, Luigi Marongiu > <marongiu.lu...@gmail.com> wrote: > > > > > > > > > > > > > > I have fitted a glm model to some data; how can I find the inverse > function of this model? Since I don't know the y=f(x) implemented by > glm (this works under the hood), I can't define a f⁻¹(y). > Is there an R function that can find the inverse of a glm model? > Thank you. > > > > > > > The working example is: > `V = c(120.64, 66.14, 34.87, 27.11, 8.87, -5.8, 4.52, -7.16, -17.39, -14.29, > -20.26, -14.99, -21.05, -20.64, -8.03, -21.56, -1.28, 15.01, 75.26, 191.76, > 455.09, 985.96, 1825.59, 2908.08, 3993.18, 5059.94, 6071.93, 6986.32, > 7796.01, 8502.25, 9111.46, 9638.01, 10077.19, 10452.02, 10751.81, 11017.49, > 11240.37, 11427.47, 11570.07, 11684.96, 11781.77, 11863.35, 11927.44, > 11980.81, 12021.88, 12058.35, 12100.63, 12133.57, 12148.89, 12137.09) df = > data.frame(Cycles = 1:35, Values = V[1:35]) M = max(df$Values) df$Norm = > df$Values/M df$Norm[df$Norm<0] = 0 b_model = glm(Norm ~ Cycles, data=df, > family=binomial) plot(Norm ~ Cycles, df, main="Normalized view", > xlab=expression(bold("Time")), ylab=expression(bold("Signal (normalized)")), > type="p", col="cyan") lines(b_model$fitted.values ~ df$Cycles, col="blue", > lwd=2)` > > > > > > > ______________________________________________ > R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide https://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Best regards, Luigi ______________________________________________ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide https://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
