# emaxnls

The **emaxnls** package provides tools for nonlinear least squares
estimation for Emax regression models. It supplies a clean interface for
specifying Emax regression models with covariates, using
[`nls()`](https://rdrr.io/r/stats/nls.html) as the underlying tool for
fitting the model. It supports least squares estimation using the
Gauss-Newton algorithm, the Levenberg-Marquardt algorithm (via
[`minpack.lm::nls.lm()`](https://rdrr.io/pkg/minpack.lm/man/nls.lm.html))
and the ‘nl2sol’ algorithm from the Port library. Continuous and binary
response variables are both supported, with an iterative reweighted
least squares method used to produce estimates in the binary case.

## Installation

You can install the latest CRAN release with:

``` R
install.packages("emaxnls")
```

Alternatively, you can install the development version of emaxnls from
[GitHub](https://github.com/) with:

``` R
# install.packages("pak")
pak::pak("djnavarro/emaxnls")
```

## Minimal example

An Emax regression model for continuous response variables is estimated
using the
[`emax_nls()`](https://emaxnls.djnavarro.net/reference/emax_nls.md)
function, using the `structural_model` argument to specify the exposure
variable and the response variable, and the `covariate_model` argument
to specify covariates to be estimated for structural parameters (e.g.,
E0, Emax).

``` R
library(tibble)
library(emaxnls)
set.seed(123)

# a synthetic data set bundled with the package
emax_df
#> # A tibble: 400 × 12
#>       id  dose  exp_1  exp_2 rsp_1 rsp_2 cnt_a cnt_b cnt_c bin_d bin_e cat_f
#>    <int> <dbl>  <dbl>  <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <fct>
#>  1     1   200 12332. 13004. 15.7      1  3.85  5.89  4.31     1     1 grp 1
#>  2     2   300 18232. 17244. 15.3      1  4.78  7.25  3.73     1     1 grp 1
#>  3     3     0     0      0   5.65     0  1.22  9.24  2.41     1     1 grp 1
#>  4     4   200  9394.  8839. 12.5      0  2.68  7.14  3.76     1     1 grp 2
#>  5     5   200  7088.  9827. 13.2      1  4.27  5.57  9.05     0     1 grp 2
#>  6     6   300 30402. 28483. 16.8      1  6.09  6.08  4.62     0     1 grp 1
#>  7     7   300 21679. 17137. 17.4      1  7.5   8.1   2.08     0     1 grp 3
#>  8     8   100 15506. 13377. 15.9      0  3.65  6.89  3.56     0     1 grp 1
#>  9     9     0     0      0   7.3      0  4.84  3.77  7.44     0     1 grp 2
#> 10    10   200  5331.  5251. 12.8      1  4.45  3.42  1.66     1     0 grp 3
#> # ℹ 390 more rows

# estimate parameters for an Emax regression with covariates
mod_c <- emax_nls(
  structural_model = rsp_1 ~ exp_1, # specify the response and exposure variables
  covariate_model = list(
    E0 ~ cnt_a,  # add a covariate on the E0 intercept parameter
    Emax ~ 1,    # no covariates on Emax
    logEC50 ~ 1  # no covariates on logEC50
  ), 
  data = emax_df
)

mod_c
#> Structural model:
#> 
#>   Exposure:       exp_1 
#>   Response:       rsp_1 
#>   Emax type:      hyperbolic 
#>   Response type:  continuous
#> 
#> Covariate model:
#> 
#>   E0:       E0 ~ cnt_a 
#>   Emax:     Emax ~ 1 
#>   logEC50:  logEC50 ~ 1 
#> 
#> Model fit:
#> 
#>   Observations:         400 
#>   Residual df:          396 
#>   Residual std. error:  0.5108 
#>   AIC:                  603.6431 
#> 
#> Coefficients (95% CI):
#> 
#>   label             estimate std_error lower  upper
#> 1 E0_cnt_a             0.486    0.0116 0.463  0.509
#> 2 E0_Intercept         5.05     0.0759 4.91   5.20 
#> 3 Emax_Intercept       9.97     0.112  9.75  10.2  
#> 4 logEC50_Intercept    8.27     0.0394 8.19   8.35 
#> 
#> Use summary() for hypothesis tests.

# hypothesis tests produced by summary()
summary(mod_c)
#> # A tibble: 4 × 7
#>   label             estimate std_error t_statistic    p_value ci_lower ci_upper
#>   <chr>                <dbl>     <dbl>       <dbl>      <dbl>    <dbl>    <dbl>
#> 1 E0_cnt_a             0.486    0.0116        42.1  3.63e-148    0.463    0.509
#> 2 E0_Intercept         5.05     0.0759        66.6  4.16e-217    4.91     5.20 
#> 3 Emax_Intercept       9.97     0.112         89.3  2.11e-264    9.75    10.2  
#> 4 logEC50_Intercept    8.27     0.0394        NA   NA            8.19     8.35
```

## Logistic Emax regression

For binary outcomes, the package provides the
[`emax_logistic()`](https://emaxnls.djnavarro.net/reference/emax_logistic.md)
function, which estimates the parameters of a logistic Emax model. Under
this model a logit-link function and binomial family is assumed, and the
parameters are estimated using an iterative reweighted nonlinear least
squares procedure.

``` R
# estimate parameters for a logistic Emax regression with covariates
mod_b <- emax_logistic(
  structural_model = rsp_2 ~ exp_1, # specify the response and exposure variables
  covariate_model = list(
    E0 ~ cnt_a,  # add a covariate on the E0 intercept parameter
    Emax ~ 1,    # no covariates on Emax
    logEC50 ~ 1  # no covariates on logEC50
  ), 
  data = emax_df
)

mod_b
#> Structural model:
#> 
#>   Exposure:       exp_1 
#>   Response:       rsp_2 
#>   Emax type:      hyperbolic 
#>   Response type:  binary (logit link)
#> 
#> Covariate model:
#> 
#>   E0:       E0 ~ cnt_a 
#>   Emax:     Emax ~ 1 
#>   logEC50:  logEC50 ~ 1 
#> 
#> Model fit:
#> 
#>   Observations:  400 
#>   Residual df:   396 
#>   Deviance:      331.4698 
#>   AIC:           339.4698 
#> 
#> Coefficients (95% CI):
#> 
#>   label             estimate std_error  lower  upper
#> 1 E0_cnt_a             0.659    0.0800  0.501  0.816
#> 2 E0_Intercept        -5.00     0.578  -6.14  -3.87 
#> 3 Emax_Intercept       8.12     2.27    5.08  17.6  
#> 4 logEC50_Intercept    9.78     0.518   8.89  11.0  
#> 
#> Use summary() for hypothesis tests.

# hypothesis tests produced by summary()
summary(mod_b)
#> # A tibble: 4 × 7
#>   label             estimate std_error z_statistic   p_value ci_lower ci_upper
#>   <chr>                <dbl>     <dbl>       <dbl>     <dbl>    <dbl>    <dbl>
#> 1 E0_cnt_a             0.659    0.0800        8.24  1.79e-16    0.501    0.816
#> 2 E0_Intercept        -5.00     0.578        -8.64  5.43e-18   -6.14    -3.87 
#> 3 Emax_Intercept       8.12     2.27          3.58  3.45e- 4    5.08    17.6  
#> 4 logEC50_Intercept    9.78     0.518        NA    NA           8.89    11.0
```

## Stepwise covariate modelling

The package supports stepwise covariate modelling via forward addition
and backward elimination. The
[`emax_scm_forward()`](https://emaxnls.djnavarro.net/reference/emax_scm.md)
function supports forward addition, the
[`emax_scm_backward()`](https://emaxnls.djnavarro.net/reference/emax_scm.md)
function supports backward elimination, and the syntax is designed to
allow forward-backward procedures by piping a base model to
[`emax_scm_forward()`](https://emaxnls.djnavarro.net/reference/emax_scm.md)
and then to
[`emax_scm_backward()`](https://emaxnls.djnavarro.net/reference/emax_scm.md).

``` R
# base model with no covariates
base_model <- emax_nls(
  structural_model = rsp_1 ~ exp_1, 
  covariate_model = list(E0 ~ 1, Emax ~ 1, logEC50 ~ 1), 
  data = emax_df
)

# list of possible consider for E0 and Emax
covariate_list <- list(
  E0 = c("cnt_a", "cnt_b", "cnt_c", "bin_d", "bin_e"),
  Emax = c("cnt_a", "cnt_b", "cnt_c", "bin_d", "bin_e")
)

# stepwise covariate modelling with a forward step and a backward step
final_mod <- base_model |> 
  emax_scm_forward(candidates = covariate_list, threshold = .01) |> 
  emax_scm_backward(candidates = covariate_list, threshold = .001) 

# extract the complete history of all models tested during the 
# stepwise covariate modelling procedure
emax_scm_history(final_mod)
#> # A tibble: 22 × 13
#>    iteration attempt step       criterion action term_tested  model_tested   model_converged convergence_reason
#>        <int>   <int> <chr>      <chr>     <chr>  <chr>        <chr>          <lgl>           <chr>             
#>  1         0       0 base model <NA>      <NA>   <NA>         E0 ~ 1, Emax … TRUE            converged         
#>  2         1       1 forward    p-value   add    E0 ~ cnt_c   E0 ~ 1 + cnt_… TRUE            converged         
#>  3         1       2 forward    p-value   add    Emax ~ bin_e E0 ~ 1, Emax … TRUE            converged         
#>  4         1       3 forward    p-value   add    E0 ~ cnt_b   E0 ~ 1 + cnt_… TRUE            converged         
#>  5         1       4 forward    p-value   add    Emax ~ cnt_c E0 ~ 1, Emax … TRUE            converged         
#>  6         1       5 forward    p-value   add    Emax ~ cnt_a E0 ~ 1, Emax … TRUE            converged         
#>  7         1       6 forward    p-value   add    Emax ~ bin_d E0 ~ 1, Emax … TRUE            converged         
#>  8         1       7 forward    p-value   add    E0 ~ cnt_a   E0 ~ 1 + cnt_… TRUE            converged         
#>  9         1       8 forward    p-value   add    Emax ~ cnt_b E0 ~ 1, Emax … TRUE            converged         
#> 10         1       9 forward    p-value   add    E0 ~ bin_e   E0 ~ 1 + bin_… TRUE            converged         
#> # ℹ 12 more rows
#> # ℹ 4 more variables: term_p_value <dbl>, model_aic <dbl>, model_bic <dbl>, model_updated <lgl>

# show the final model
final_mod
#> Structural model:
#> 
#>   Exposure:       exp_1 
#>   Response:       rsp_1 
#>   Emax type:      hyperbolic 
#>   Response type:  continuous
#> 
#> Covariate model:
#> 
#>   E0:       E0 ~ 1 + cnt_a 
#>   Emax:     Emax ~ 1 
#>   logEC50:  logEC50 ~ 1 
#> 
#> Model fit:
#> 
#>   Observations:         400 
#>   Residual df:          396 
#>   Residual std. error:  0.5108 
#>   AIC:                  603.6431 
#> 
#> Coefficients (95% CI):
#> 
#>   label             estimate std_error lower  upper
#> 1 E0_cnt_a             0.486    0.0116 0.463  0.509
#> 2 E0_Intercept         5.05     0.0759 4.91   5.20 
#> 3 Emax_Intercept       9.97     0.112  9.75  10.2  
#> 4 logEC50_Intercept    8.27     0.0394 8.19   8.35 
#> 
#> Use summary() for hypothesis tests.
```

## Simulation

The package also provides tools to assist in model-based simulations,
using the
[`simulate()`](https://emaxnls.djnavarro.net/reference/simulate.md)
function. A simple example is shown below. See the function
documentation for more details.

``` R
simulate(final_mod, nsim = 1)
#> # A tibble: 400 × 11
#>    dat_id sim_id    mu   val E0_cnt_a E0_Intercept Emax_Intercept logEC50_Intercept rsp_1  exp_1 cnt_a
#>     <int>  <int> <dbl> <dbl>    <dbl>        <dbl>          <dbl>             <dbl> <dbl>  <dbl> <dbl>
#>  1      1      1 14.5  13.9     0.486         5.01           9.95              8.25 15.7  12332.  3.85
#>  2      2      1 15.6  15.5     0.486         5.01           9.95              8.25 15.3  18232.  4.78
#>  3      3      1  5.60  6.15    0.486         5.01           9.95              8.25  5.65     0   1.22
#>  4      4      1 13.4  13.3     0.486         5.01           9.95              8.25 12.5   9394.  2.68
#>  5      5      1 13.6  13.0     0.486         5.01           9.95              8.25 13.2   7088.  4.27
#>  6      6      1 16.8  16.4     0.486         5.01           9.95              8.25 16.8  30402.  6.09
#>  7      7      1 17.1  17.5     0.486         5.01           9.95              8.25 17.4  21679.  7.5 
#>  8      8      1 14.8  14.6     0.486         5.01           9.95              8.25 15.9  15506.  3.65
#>  9      9      1  7.36  6.69    0.486         5.01           9.95              8.25  7.3      0   4.84
#> 10     10      1 13.0  12.7     0.486         5.01           9.95              8.25 12.8   5331.  4.45
#> # ℹ 390 more rows
```

# Package index

## Emax regression

Build continuous response Emax regression models

- [`emax_nls()`](https://emaxnls.djnavarro.net/reference/emax_nls.md) :
  Estimate parameters for an Emax regression model
- [`emax_nls_options()`](https://emaxnls.djnavarro.net/reference/emax_nls_options.md)
  : Settings used to estimate Emax model
- [`emax_nls_init()`](https://emaxnls.djnavarro.net/reference/emax_nls_init.md)
  : Construct an initial guess for the Emax model parameters

## Emax logistic regression

Build binary response Emax logistic regression models

- [`emax_logistic()`](https://emaxnls.djnavarro.net/reference/emax_logistic.md)
  : Estimate parameters for a logistic Emax regression model
- [`emax_logistic_options()`](https://emaxnls.djnavarro.net/reference/emax_logistic_options.md)
  : Settings used to estimate a logistic Emax model
- [`emax_logistic_init()`](https://emaxnls.djnavarro.net/reference/emax_logistic_init.md)
  : Construct an initial guess for logistic Emax model parameters

## Covariate selection

Stepwise covariate modeling for Emax regression models

- [`emax_scm_forward()`](https://emaxnls.djnavarro.net/reference/emax_scm.md)
  [`emax_scm_backward()`](https://emaxnls.djnavarro.net/reference/emax_scm.md)
  [`emax_scm_history()`](https://emaxnls.djnavarro.net/reference/emax_scm.md)
  : Stepwise covariate modelling for Emax regression

## Methods for Emax regression objects

S3 methods for `emaxnls` and `emaxlogistic` objects

- [`print(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/print.md)
  [`print(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/print.md)
  : Print an Emax regression model object
- [`summary(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/summary.md)
  [`summary(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/summary.md)
  : Summary of an Emax regression model
- [`coef(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/coef.emaxnls.md)
  : Coefficients for an Emax regression
- [`vcov(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/vcov.emaxnls.md)
  : Variance-covariance matrix for an Emax regression
- [`confint(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/confint.emaxnls.md)
  : Confidence intervals for Emax regression model parameters
- [`nobs(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/nobs.emaxnls.md)
  : Number of observations for an Emax regression model
- [`sigma(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/sigma.emaxnls.md)
  : Residual standard deviation for an Emax regression model
- [`residuals(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/residuals.md)
  [`residuals(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/residuals.md)
  : Residuals for an Emax regression model
- [`fitted(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/fitted.md)
  [`fitted(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/fitted.md)
  : Fitted values for an Emax regression model
- [`predict(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/predict.md)
  [`predict(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/predict.md)
  : Predicting from Emax regression models
- [`simulate(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/simulate.md)
  [`simulate(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/simulate.md)
  : Simulate responses from an Emax regression model
- [`logLik(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/logLik.md)
  [`logLik(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/logLik.md)
  : Log-likelihood for an Emax regression model
- [`AIC(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/AIC.md)
  [`BIC(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/AIC.md)
  [`AIC(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/AIC.md)
  [`BIC(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/AIC.md)
  : Akaike information criterion / Bayesian information criterion
- [`anova(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/anova.md)
  [`anova(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/anova.md)
  : Analysis of variance for Emax regression models
- [`deviance(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/deviance.md)
  [`deviance(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/deviance.md)
  : Model deviance for an Emax regression model
- [`df.residual(`*`<emaxlogistic>`*`)`](https://emaxnls.djnavarro.net/reference/df.residual.md)
  [`df.residual(`*`<emaxnls>`*`)`](https://emaxnls.djnavarro.net/reference/df.residual.md)
  : Residual degrees of freedom for an Emax regression model

## Other

Other exported functions and objects

- [`emax_df`](https://emaxnls.djnavarro.net/reference/emax_df.md) :
  Sample simulated data for Emax exposure-response models with
  covariates.
- [`emax_converged()`](https://emaxnls.djnavarro.net/reference/emax_converged.md)
  : Check Emax regression model for convergence status
- [`emax_add_term()`](https://emaxnls.djnavarro.net/reference/emax_update.md)
  [`emax_remove_term()`](https://emaxnls.djnavarro.net/reference/emax_update.md)
  : Add or remove a covariate term from an Emax regression
- [`emax_fun()`](https://emaxnls.djnavarro.net/reference/emax_fun.md) :
  Construct Emax prediction function from model object

# Articles

### Articles

- [Fitting Emax models for continuous
  outcomes](https://emaxnls.djnavarro.net/articles/fitting-emax-models.md):
- [Fitting Emax models for binary
  outcomes](https://emaxnls.djnavarro.net/articles/fitting-logistic-emax-models.md):
- [Simulating from Emax
  models](https://emaxnls.djnavarro.net/articles/simulating-from-emax-models.md):
- [Stepwise covariate
  modelling](https://emaxnls.djnavarro.net/articles/stepwise-covariate-modelling.md):
- [Visualising exposure-response models with
  erplots](https://emaxnls.djnavarro.net/articles/erplots-integration.md):