All of your prediction interval questions about logistic regression are answered in the first attached screenshot search. The second attachment tells you how to build these prediction models. Best wishes, David Booth

I don't think there is a package providing a general solution, because the solution can be arbitrarily complex (high-dimensional sets, depending on the structural part of the model that may include any number of predictors, interactions and non-linear transformations).

For the univariate case with linear (or at least monotone) relationship between predictor and log odds the solution is quite simple. You could use optim(), optimise() or nlm() to find the x-coordinate where the lower and uppper limit of the confidence band have a given value (0.5 for instance).

But one should never underestimate the creativity and ability of R developers.... so there might in fact be some solution available for more complex models. I'm placing myself in the chain waiting for a good hint coming up :)

The attached screenshot should be of assistance in using R packages in your research. Apologies for not adding this earlier.. Best wishes to all, David Booth

Thanks David, Danke Jochen.

My collaborators and I are getting closer to a solution with these two recent papers that I've just found and which are freely available (Open Access):

1) Fisher, S. M., and M. A. Lewis. 2021. A robust and efficient algorithm to find profile likelihood confidence intervals. Statistics and Computing:31:38. The authors provide their Python codes to estimate the profile likelihood CI of a logistic regression model (both for its parameters and more interestingly its prediction).

2) Ren, X., and J. Xia. 2019. An algorithm for computing profile likelihood based pointwise confidence intervals for nonlinear dose-response models. PLoS ONE 14:e0210953. The authors provides SAS codes to estimate the profile likelihood CI of the LD50 which is analogous to the L50 in fish and is probably our best guess to start with. Trying to apply the profile likelihood to our parameter of interest (L50) seems like a good idea, as it outperformed the Wald- and Delta-based methods in terms of coverage.

Unfortunately, none of these two papers is based on R and the codes are not straightforward to interpret/transfer to R for a biologist like me, but we have something here to work with and my statistician collaborator should be able to get things work out in the end.

Looking forward still to see if a R package may exist somewhere for profile likelihood CIs of a logistic regression model prediction or if one may be released in a near future.

Why not call the python code from R?

https://cran.r-project.org/web/packages/reticulate/vignettes/calling_python.html

It might be a bit inefficient for simulation work, but worth a try.

Thanks Thom, we will certainly have a look at this option. Cheers!

My initial inclination would be to derive the confidence intervals by bootstrap. As long as you can write the function for the model and the subsequent L50, you can just pass this function to the bootci() function in the boot package. I don't know if this is an ideal solution, but it's a simple one.

Hi Salvatore,

Thanks for your suggestion.

My collaborator and I will be comparing in fact 7 different methods to assess the uncertainty (95% CI) in a logistic regression model prediction, i.e. the L50 in our case (i.e., the length at which a fish has a 50% chance of exhibiting mature gonads):

1- Delta method

2- Wald-based method

3- Profile likelihood (the one for which an R package would be great)

4- Non-parametric bootstrapping

5- Parametric bootstrapping (less common)

6- Fieller (1944) analytical method

7- Bayesian credible intervals

The profile likelihood function is of interest to us, as it better models the sometimes non-normal variance in the linear predictors on the link scale. Moreover, as pointed out by Royston (2007) regarding the advantage of the profile likelihood over the Wald and bootstrap approaches:

"Both examples indicate that using profile likelihood improves on normal-based CIs. The main reason is that the likelihood-ratio statistic tends to approach its asymptotic distribution more rapidly than the equivalent Wald statistic [...] In principle, one could do so also by using the bootstrap. However, for CI calculations not assuming normality, the bootstrap can become computationally expensive. The bootstrap is usually assumed valid, even in small samples, but that may not be so; the bootstrap gives only asymptotically correct results."

And Zheng et al. (2012) concluded that the profile likelihood function should be preferred over the non-parametric bootstrapping approach, as the latter performed more poorly in controlling type 1 error.

Looking forward to see what our study will reveal, but we'll need to figure out how to adequately quantify the profile likelihood CI for this.

Thanks again for your suggestion.

Royston, P. 2007. Profile likelihood for estimation and confidence intervals. The Stata Journal 7:376-387.

Zheng, H., Z. Liao, S. Liu, W. Liang, X. Zhang, and C. Ou. 2012. Comparison of three methods in estimating confidence intervals and hypothesis testing of logistic regression coefficients. Journal of Mathematical Medicine 4:393-396.

My collaborator has been able to code the profile likelihood function in R (as well as other methods) for our study of the L50 in fish and the assessment of its uncertainty. If our work get published, we will make it available (if allowed to).

Thanks to all for the input.

Cheers!

Two years ago, we published an algorithm to compute profile likelihood confidence intervals based on the likelihood function. The paper has been mentioned in the discussion already: Article A robust and efficient algorithm to find profile likelihood ...

Now I have provided a way to apply this algorithm from R (in addition to the original Python implementation). Please have a look at the updated documentation page: https://vemomoto.github.io/ci_rvm.html#usage-from-r

It is not trivial to test the package from the perspective of an end user, so if you encounter any problems, please feel free to contact me or use the bug tracker on the project's GitHub page: https://github.com/vemomoto/vemomoto/tree/master/ci_rvm

Guten Tag Samuel Matthias Fischer,

Thanks for having made your [trust-region based] profile-likelihood method available for R users. I'm confindent it's a clear improvement over the original method of Venzon & Moolgavkar (1988).

Our study about the different approaches to estimate the uncertainty around the L50 in fish is being reviewed as we speak, so unfortunately we won't be able to include it as well in our work. Know however that the original profile-likelihood did pretty well (ranked 4th out of 13 considered methods) in our comparison tests using simulations with different effect and sample sizes; much better actually than the different forms of often-used bootstrapping approaches that are more prone to accumulate numerical errors and thus generates sometimes biased CI bounds.

But to come back to the main point, we'll have a look at your method in in a near future hopefully and let you know if we've been able to use it. For instance, to compare it with the top-ranking method we have identified based on our L50 simulations derived from > 23,000 walleye females sampled over two decades.

Danke sehr!