Solved – How to find the inflection point in logistic regression using R

Question

I have conducted a binary logistic regression (with 7 predictor variables) in R using "year" as a random effect, i.e. a generalized linear mixed model (GLMM). For a few of the predictor variables, I would now like to find the inflection point, i.e. the value of the predictor variable for which the probability of the response variable being "1" equals 0.5. I would like both to produce a curve of this and to calculate the exact value, so that I won't have to read it from the curve (with the inaccuracy this would cause). A post from November 2013 made by Francis and edited by Scortchi touches on this topic (with a different purpose, though) and shares a specific code for producing such a curve (Logistic Regression and Inflection Point). However, I don't understand the 4th line of the code (newdata) and can't make the code work for my data. Is the 4th line something that has to be done in order to produce the curve? Might someone help me by uploading a more general code, or gently translating the referred one, for me to understand exactly how to fit my own data into the code? Also, if someone knows how to make R calculate the exact value, I would very much like to know. As you might be able to tell, I'm still learning R. I hope that someone can help.

Answer 1

When the logit of the probability $\pi$ of the response's being 1 is linear in the $i$th predictor $x_i$, you can write the model like this

$$\log \frac{\pi}{1-\pi} = \beta_i x_i + \beta_0 +\sum_{j\neq i} \beta_j x_j$$

where the $\beta$s are the coefficients, & the $x_j$s the other predictors. The inflection point $x_i^*$ is where

$$\left.\frac{\partial^2 \pi}{(\partial x_i)^2}\right|_{x_i=x_i^*}=0$$

When you work that out it gives

$$x_i^*=\frac{-(\beta_0 +\sum_{j\neq i} \beta_j x_j)}{\beta_i}$$

& the corresponding probability, $\pi_i^*$, doesn't depend on the coefficients or other predictor values:

$$\pi_i^*=\frac{1}{2}$$

To extract coefficient estimates in R use coefficients(my.model). A convenient way to plot the curve is to make a data frame my.curve with $x_i$ varying over a range of interest & the $x_j$s constant (this is the line of code you queried) & then use predict(mymodel, newdata=my.curve, type="response") to obtain the predicted probabilities.