delicatessen.estimating_equations.pharmacokinetics.ee_loglogistic

ee_loglogistic(theta, dose, response, robust=None, k=None)

Estimating equations for the 4 parameter log-logistic dose-response model. The log-logistic model describes the dose-response relationship in terms of four parameters: the zero-dose response, the maximum response (E-max), the dose producing half maximal effect (ED50), and steepness of the dose-response curve. The assumed model is

\[R_i = \theta_{z} + \frac{\theta_{m} - \theta_{z}}{1 + \exp\left[ \theta_{s} (\log(D_i) - \log(\theta_{50})) \right]}\]

where \(R\) is the response and \(D\) is the dose. Here, \(\theta_{z}\) is the zero-dose response, \(\theta_{m}\) is the maximum response, \(\theta_{50}\) is the dose with 50% of maximal response, and \(\theta_{s}\) is the slope. The corresponding estimating equations for this model are

\[\begin{split}\sum_{i=1}^n \begin{bmatrix} -2 (Y_i - \hat{Y}_i) (1 - 1/(1 + \rho)) \\ 2 (Y_i - \hat{Y}_i) (1 / (1 + \rho))) \\ 2 (Y_i - \hat{Y}_i) (\theta_m - \theta_z) \frac{\theta_2}{\theta_1} \frac{\rho}{(1 + \rho)^2} \\ 2 (Y_i - \hat{Y}_i) (\theta_m - \theta_z) \log(D_i / \theta_1) \frac{\rho}{(1 + \rho)^2} \end{bmatrix} = 0\end{split}\]

where \(R_i\) is the response of individual \(i\), \(D_i\) is the dose, \(\rho = \frac{D_i}{\theta_{50}}^{\theta_s}\), and \(\hat{Y_i} = \theta_z + \frac{\theta_m - \theta_z}{1+\rho}\).

Here, theta is a 1-by-4 array. The first theta corresponds to lower limit (\(\theta_z\)), the second corresponds to the upper limit (\(\theta_m\)), the third corresponds to the effective dose (ED50) (\(\theta_{50}\)), and the fourth corresponds to the steepness (\(\theta_s\)).

Note

This implementation supports models where dose increases response and dose decreases response. Depending on the relationship observed in the data, set the starting values for the lower and upper responses according (i.e., increasing has lower < upper and descreasing has `` lower > upper``). The steepness parameter should be positive for increasing and negative for decreasing. If starting parameters are not chosen well, the model may not converge or converge to nonsensical values.

Parameters

  • theta (ndarray, list, vector) – Theta in this case consists of 4 values. Outside of steepness (i.e., theta[2]), non-negative starting values are recommended for the log-logistic model.

  • dose (ndarray, list, vector) – 1-dimensional vector of n dose values.

  • response (ndarray, list, vector) – 1-dimensional vector of n response values.

  • robust (None, str, optional) – Robust loss function to control the influence of outliers. Default is None, which uses the standard (non-outlier) robust version of the model. For a complete set of options for the outlier-robust functions, see the robust_loss_functions documentation.

  • k (None, int, float, optional) – Tuning or hyperparameter for the chosen outlier-robust function. Notice that the choice of hyperparameter should depend on the chosen loss function.

Returns

Returns a 4-by-n NumPy array evaluated for the input theta.

Return type

array

Examples

Construction of a estimating equation(s) with ee_emax_sigmoid should be done similar to the following

>>> from delicatessen import MEstimator
>>> from delicatessen.data import load_inderjit
>>> from delicatessen.estimating_equations import ee_loglogistic

For demonstration, we use dose-response data from Inderjit et al. (2002), which can be loaded from delicatessen directly.

>>> d = load_inderjit()   # Loading array of data
>>> resp_data = d[:, 0]   # Response data
>>> dose_data = d[:, 1]   # Dose data

Defining psi, or the stacked estimating equations

>>> def psi(theta):
>>>     return ee_loglogistic(theta=theta, dose=dose_data, response=resp_data)

The sigmoid E-max model is harder to solve compared to other estimating equations. Namely, the root-finder is not aware of implicit bounds on the parameters. To reduce non-convergence issues, we can give the root-finder good starting values.

Here, we use some general starting values that should perform well in many cases. For the lower-bound, give the minimum response value as the initial. For ED50, give the mid-point between the maximum response and the minimum response. From the scatter plot of the dose-response, we see that the response increases as dose increases. So, the starting value for steepness should be positive. For the upper-bound, give the maximum response value as the initial.

>>> estr = MEstimator(psi, init=[np.min(resp_data),
>>>                              np.max(resp_data),
>>>                              (np.max(resp_data)+np.min(resp_data)) / 2,
>>>                              (np.max(resp_data)+np.min(resp_data)) / 2])
>>> estr.estimate(solver='lm')

Inspecting the parameter estimates, variance, and confidence intervals

>>> estr.theta
>>> estr.variance
>>> estr.confidence_intervals()

Inspecting the parameter estimates

>>> estr.theta[0]    # lower limit
>>> estr.theta[1]    # upper limit
>>> estr.theta[2]    # ED50
>>> estr.theta[3]    # steepness

One can also use ee_emax_sigmoid to estimate 3-parameter and 2-parameter log-logistic models. For these models, a constant is added to the input array and only a subset of the output estimating equations are used. The following is an example of how to estimate a 3 parameter log-logistic model, which assumes that the lower limit of the response is zero (this makes sense in the context of this application).

>>> def psi(theta):
>>>     theta = [0, ] + list(theta)    # Setting lower-limit manually to zero
>>>     ee = ee_loglogistic(theta=theta, dose=dose_data, response=resp_data)
>>>     return ee[1:, :]               # Return all estimating equations besides lower limit

Since we are no longer root-finding for the lower limit (i.e., it is manually set to zero in psi), only three starting values need to be provided

>>> estr = MEstimator(psi, init=[(np.max(resp_data)+np.min(resp_data)) / 2,
>>>                              (np.max(resp_data)+np.min(resp_data)) / 2,
>>>                              np.max(resp_data)])
>>> estr.estimate(solver='lm')

Inspecting the parameter estimates, variance, and confidence intervals

>>> estr.theta
>>> estr.variance
>>> estr.confidence_intervals()

Inspecting the parameter estimates (notice that the indices differ from the usual sigmoid E-max, since we dropped the lower-limit estimating equation).

>>> estr.theta[0]    # upper limit
>>> estr.theta[1]    # ED50
>>> estr.theta[2]    # steepness

References

An H, Justin TL, Aubrey GB, Marron JS, & Dittmer DP. (2019). dr4pl: A Stable Convergence Algorithm for the 4 Parameter Logistic Model. R J., 11(2), 171.

Fomenko I, Durst M, & Balaban D. (2006). Robust regression for high throughput drug screening. Computer Methods and Programs in Biomedicine, 82(1), 31-37.

Inderjit, Streibig JC, & Olofsdotter M. (2002). Joint action of phenolic acid mixtures and its significance in allelopathy research. Physiologia Plantarum, 114(3), 422-428.

Lim C, Sen PK, & Peddada SD. (2013). Robust nonlinear regression in applications. Journal of the Indian Society of Agricultural Statistics. Indian Society of Agricultural Statistics, 67(2), 215.

Ritz C, Baty F, Streibig JC, & Gerhard D. (2015). Dose-response analysis using R. PloS One, 10(12), e0146021.