Firth's Bias-Reduced Logistic Regression (logistf)

Usage

strategy_logistf(
  modelstr,
  model_name = "logistf",
  report_columns = c("statistic", "p.value", "p.value.adjusted", "moderated.p.value",
    "moderated.p.value.adjusted"),
  test = "Chisq"
)

strategy_lmer(
  modelstr,
  model_name = "Model",
  report_columns = c("statistic", "p.value", "p.value.adjusted", "moderated.p.value",
    "moderated.p.value.adjusted")
)

strategy_lm(
  modelstr,
  model_name = "Model",
  report_columns = c("statistic", "p.value", "p.value.adjusted", "moderated.p.value",
    "moderated.p.value.adjusted")
)

strategy_rlm(
  modelstr,
  model_name = "Model",
  report_columns = c("statistic", "p.value", "p.value.adjusted", "moderated.p.value",
    "moderated.p.value.adjusted")
)

strategy_glm(
  modelstr,
  model_name = "Model",
  test = "Chisq",
  family = stats::binomial,
  multiplier = 1,
  offset = 1,
  report_columns = c("statistic", "p.value", "p.value.adjusted", "moderated.p.value",
    "moderated.p.value.adjusted")
)

Arguments

modelstr: model formula
model_name: name of model
report_columns: columns to report
test: type of test statistic passed to anova (e.g. "Chisq")
family: either binomial or quasibinomial
multiplier: for tuning default is 1.
offset: offset added to contingency table cells to avoid perfect separation

Value

list with model function, contrast computation function etc.

Other modelling: Contrasts, ContrastsFirth, ContrastsLimma, ContrastsMissing, ContrastsModerated, ContrastsModeratedDEqMS, ContrastsPlotter, ContrastsProDA, ContrastsROPECA, ContrastsTable, INTERNAL_FUNCTIONS_BY_FAMILY, LR_test(), Model, ModelFirth, ModelLimma, build_model(), build_model_limma(), build_model_logistf(), contrasts_fisher_exact(), get_anova_df(), get_complete_model_fit(), get_p_values_pbeta(), group_label(), isSingular_lm(), linfct_all_possible_contrasts(), linfct_factors_contrasts(), linfct_from_model(), linfct_matrix_contrasts(), merge_contrasts_results(), model_analyse(), model_summary(), moderated_p_deqms(), moderated_p_deqms_long(), moderated_p_limma(), moderated_p_limma_long(), my_contest(), my_contrast(), my_contrast_V1(), my_contrast_V2(), my_glht(), pivot_model_contrasts_2_Wide(), plot_lmer_peptide_predictions(), sim_build_models_lm(), sim_build_models_lmer(), sim_build_models_logistf(), sim_make_model_lm(), sim_make_model_lmer(), strategy_limma(), summary_ROPECA_median_p.scaled()

Examples


tmp <- strategy_logistf("bin_resp ~ condition", model_name = "parallel design")
tmp$model_fun(get_formula = TRUE)
#> bin_resp ~ condition
#> <environment: 0x562abf0d7130>
tmp$isSingular
#> function (m) 
#> {
#>     anyNA <- any(is.na(coefficients(m)))
#>     if (anyNA) {
#>         return(TRUE)
#>     }
#>     else {
#>         if (df_residual_logistf(m) >= 2) {
#>             return(FALSE)
#>         }
#>         return(TRUE)
#>     }
#> }
#> <bytecode: 0x562ab4d8ebb8>
#> <environment: 0x562abf0d7130>

istar <- prolfqua::sim_lfq_data_peptide_config(Nprot = 10, with_missing = TRUE,
  weight_missing = 0.5, seed = 3)
#> creating sampleName from fileName column
#> completing cases
#> completing cases done
#> setup done
istar$data <- encode_bin_resp(istar$data, istar$config)
#> completing cases
istar <- LFQData$new(istar$data, istar$config)
df <- istar$summarize_hierarchy()
df2 <- df[df[[ncol(df)]] > 1,  ]
istar2 <- istar$get_subset(df2)
#> Joining with `by = join_by(protein_Id)`
istar2$data |>
dplyr::group_by(protein_Id) |>
 tidyr::nest() -> nestProtein
modelFunction <- strategy_logistf("bin_resp ~ group_ + peptide_Id", model_name = "random_example")
modelFunction$model_fun(nestProtein$data[[1]])
#> logistf::logistf(formula = formula, data = DFT, weights = Freq)
#> Model fitted by Penalized ML
#> Confidence intervals and p-values by Profile Likelihood 
#> 
#> Coefficients:
#>        (Intercept)            group_B         group_Ctrl peptide_IdFLq7LKTq 
#>       2.101899e+00       6.773389e-01      -6.663876e-01       3.440276e-16 
#> peptide_IdJYhOpuPH peptide_IdLiw5EMKP peptide_IdVcatZJTa peptide_IdjrLUqOjg 
#>      -1.068335e+00       3.432936e-16      -1.068335e+00       1.197563e+00 
#> peptide_Idq2jTaC1y 
#>      -1.445323e+00 
#> 
#> Likelihood ratio test=10.07111 on 8 df, p=0.2600707, n=84
#> 
modelFunction$model_fun(nestProtein$data[[4]])
#> logistf::logistf(formula = formula, data = DFT, weights = Freq)
#> Model fitted by Penalized ML
#> Confidence intervals and p-values by Profile Likelihood 
#> 
#> Coefficients:
#>        (Intercept)            group_B         group_Ctrl peptide_IdWcAw5ozd 
#>       8.375544e-03       8.360400e-11       9.365126e-01      -3.147690e-01 
#> peptide_IdgdnXrza3 peptide_IdxvlVt88v 
#>      -2.035802e-09      -6.303721e-01 
#> 
#> Likelihood ratio test=3.276542 on 5 df, p=0.657435, n=48
#> 


istar <- prolfqua::sim_lfq_data_peptide_config(Nprot = 10, with_missing = FALSE)
#> creating sampleName from fileName column
#> completing cases
#> completing cases done
#> setup done
istar <- prolfqua::LFQData$new(istar$data,istar$config)
istar$data <- istar$data |> dplyr::group_by(protein_Id) |>
dplyr::mutate(abundanceC = abundance - mean(abundance)) |> dplyr::ungroup()
istar$factors()
#> # A tibble: 12 × 3
#>    sample  sampleName group_
#>    <chr>   <chr>      <chr> 
#>  1 A_V1    A_V1       A     
#>  2 A_V2    A_V2       A     
#>  3 A_V3    A_V3       A     
#>  4 A_V4    A_V4       A     
#>  5 B_V1    B_V1       B     
#>  6 B_V2    B_V2       B     
#>  7 B_V3    B_V3       B     
#>  8 B_V4    B_V4       B     
#>  9 Ctrl_V1 Ctrl_V1    Ctrl  
#> 10 Ctrl_V2 Ctrl_V2    Ctrl  
#> 11 Ctrl_V3 Ctrl_V3    Ctrl  
#> 12 Ctrl_V4 Ctrl_V4    Ctrl  
modelFunction <- strategy_lmer("abundanceC ~ group_ + (1|peptide_Id) ",
  model_name = "random_example")
mod <- build_model(
 istar,
 modelFunction)
#> Warning: There were 4 warnings in `dplyr::mutate()`.
#> The first warning was:
#> ℹ In argument: `linear_model = purrr::map(data, model_strategy$model_fun, pb =
#>   pb)`.
#> ℹ In group 2: `protein_Id = "7cbcrd~5725"`.
#> Caused by warning in `value[[3L]]()`:
#> ! WARN :Error: grouping factors must have > 1 sampled level
#> ℹ Run `dplyr::last_dplyr_warnings()` to see the 3 remaining warnings.
#> Joining with `by = join_by(protein_Id)`
sum(mod$modelDF$exists_lmer)
#> [1] 6
sum(mod$modelDF$isSingular, na.rm=TRUE)
#> [1] 0



tmp <- strategy_lm("Intensity ~ condition", model_name = "parallel design")
tmp$model_fun(get_formula = TRUE)
#> Intensity ~ condition
#> <environment: 0x562ac845a510>
tmp$isSingular
#> function (m) 
#> {
#>     anyNA <- any(is.na(coefficients(m)))
#>     if (anyNA) {
#>         return(TRUE)
#>     }
#>     else {
#>         if (df.residual(m) >= 2) {
#>             return(FALSE)
#>         }
#>         return(TRUE)
#>     }
#> }
#> <bytecode: 0x562aae5b1db8>
#> <environment: namespace:prolfqua>
tmp <- strategy_rlm("Intensity ~ condition", model_name = "parallel design")
tmp$model_fun(get_formula = TRUE)
#> Intensity ~ condition
#> <environment: 0x562ab894c698>
tmp$isSingular
#> function (m) 
#> {
#>     anyNA <- any(is.na(coefficients(m)))
#>     if (anyNA) {
#>         return(TRUE)
#>     }
#>     else {
#>         if (df.residual(m) >= 2) {
#>             return(FALSE)
#>         }
#>         return(TRUE)
#>     }
#> }
#> <bytecode: 0x562aae5b1db8>
#> <environment: namespace:prolfqua>
tmp <- strategy_glm("Intensity ~ condition", model_name = "parallel design")
tmp$model_fun(get_formula = TRUE)
#> Intensity ~ condition
#> <environment: 0x562ac1a19130>
tmp$isSingular
#> function (m) 
#> {
#>     anyNA <- any(is.na(coefficients(m)))
#>     if (anyNA) {
#>         return(TRUE)
#>     }
#>     else {
#>         if (df.residual(m) >= 2) {
#>             return(FALSE)
#>         }
#>         return(TRUE)
#>     }
#> }
#> <bytecode: 0x562aae5b1db8>
#> <environment: namespace:prolfqua>

Usage

Arguments

Value

See also

Examples