Skip to content

Commit

Permalink
draft p_direction
Browse files Browse the repository at this point in the history
  • Loading branch information
strengejacke committed Sep 4, 2024
1 parent b373476 commit 0d60d99
Show file tree
Hide file tree
Showing 15 changed files with 563 additions and 33 deletions.
2 changes: 1 addition & 1 deletion DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
Type: Package
Package: parameters
Title: Processing of Model Parameters
Version: 0.22.2.2
Version: 0.22.2.3
Authors@R:
c(person(given = "Daniel",
family = "Lüdecke",
Expand Down
19 changes: 19 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -350,6 +350,22 @@ S3method(model_parameters,zeroinfl)
S3method(model_parameters,zoo)
S3method(p_calibrate,default)
S3method(p_calibrate,numeric)
S3method(p_direction,coxph)
S3method(p_direction,feis)
S3method(p_direction,felm)
S3method(p_direction,gee)
S3method(p_direction,glm)
S3method(p_direction,glmmTMB)
S3method(p_direction,gls)
S3method(p_direction,hurdle)
S3method(p_direction,lm)
S3method(p_direction,lme)
S3method(p_direction,merMod)
S3method(p_direction,mixed)
S3method(p_direction,rma)
S3method(p_direction,svyglm)
S3method(p_direction,wbm)
S3method(p_direction,zeroinfl)
S3method(p_significance,coxph)
S3method(p_significance,feis)
S3method(p_significance,felm)
Expand Down Expand Up @@ -563,6 +579,7 @@ S3method(print,n_clusters_hclust)
S3method(print,n_clusters_silhouette)
S3method(print,n_factors)
S3method(print,p_calibrate)
S3method(print,p_direction_lm)
S3method(print,p_significance_lm)
S3method(print,parameters_brms_meta)
S3method(print,parameters_clusters)
Expand Down Expand Up @@ -930,6 +947,7 @@ export(n_components)
export(n_factors)
export(n_parameters)
export(p_calibrate)
export(p_direction)
export(p_function)
export(p_significance)
export(p_value)
Expand Down Expand Up @@ -969,6 +987,7 @@ export(supported_models)
export(visualisation_recipe)
importFrom(bayestestR,ci)
importFrom(bayestestR,equivalence_test)
importFrom(bayestestR,p_direction)
importFrom(bayestestR,p_significance)
importFrom(datawizard,demean)
importFrom(datawizard,describe_distribution)
Expand Down
4 changes: 4 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,10 @@
* Used more accurate analytic approach to calculate normal distributions for
the SGPV in `equivalence_test()` and used in `p_significance()`.

* Added `p_direction()` methods for frequentist models. This is a convenient
way to test the direction of the effect, which formlery was already (and still
is) possible with `pd = TRUE` in `model_parameters()`.

# parameters 0.22.2

## New supported models
Expand Down
4 changes: 3 additions & 1 deletion R/equivalence_test.R
Original file line number Diff line number Diff line change
Expand Up @@ -105,7 +105,9 @@ bayestestR::equivalence_test
#' represents _the proportion of data-supported hypotheses that are also null
#' hypotheses_ _(Blume et al. 2018, Lakens and Delacre 2020)_. It represents the
#' proportion of the _full_ confidence interval range (assuming a normally
#' distributed, equal-tailed interval) that is inside the ROPE.
#' distributed, equal-tailed interval) that is inside the ROPE. The SGPV ranges
#' from zero to one. Higher values indicate that the effect is more likely to be
#' practically equivalent ("not of interest").
#'
#' Note that the assumed interval, which is used to calculate the SGPV, is an
#' estimation of the _full interval_ based on the chosen confidence level. For
Expand Down
198 changes: 198 additions & 0 deletions R/p_direction.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,198 @@
#' @importFrom bayestestR p_direction
#' @export
bayestestR::p_direction


#' @title Probability of Direction (pd)
#' @name p_direction.lm
#'
#' @description Compute the **Probability of Direction** (*pd*, also known as
#' the Maximum Probability of Effect - *MPE*). This can be interpreted as the
#' probability that a parameter (described by its full confidence, or
#' "compatibility" interval) is strictly positive or negative (whichever is the
#' most probable). Although differently expressed, this index is fairly similar
#' (i.e., is strongly correlated) to the frequentist *p-value* (see 'Details').
#'
#' @param x A statistical model.
#' @inheritParams bayestestR::p_direction
#' @inheritParams model_parameters.default
#' @param ... Arguments passed to other methods, e.g. `ci()`.
#'
#' @seealso See also [`equivalence_test()`], [`p_function()`] and
#' [`p_significance()`] for functions related to checking effect existence and
#' significance.
#'
#' @inheritSection bayestestR::p_direction What is the *pd*?
#'
#' @inheritSection bayestestR::p_direction Relationship with the p-value
#'
#' @inheritSection bayestestR::p_direction Possible Range of Values
#'
#' @inheritSection model_parameters Statistical inference - how to quantify evidence
#'
#' @references
#'
#' - Amrhein, V., Korner-Nievergelt, F., and Roth, T. (2017). The earth is
#' flat (p > 0.05): Significance thresholds and the crisis of unreplicable
#' research. PeerJ, 5, e3544. \doi{10.7717/peerj.3544}
#'
#' - Greenland S, Rafi Z, Matthews R, Higgs M. To Aid Scientific Inference,
#' Emphasize Unconditional Compatibility Descriptions of Statistics. (2022)
#' https://arxiv.org/abs/1909.08583v7 (Accessed November 10, 2022)
#'
#' - Lakens, D. (2024). Improving Your Statistical Inferences (Version v1.5.1).
#' Retrieved from https://lakens.github.io/statistical_inferences/.
#' \doi{10.5281/ZENODO.6409077}
#'
#' - Lakens, D., Scheel, A. M., and Isager, P. M. (2018). Equivalence Testing
#' for Psychological Research: A Tutorial. Advances in Methods and Practices
#' in Psychological Science, 1(2), 259–269. \doi{10.1177/2515245918770963}
#'
#' - Makowski, D., Ben-Shachar, M. S., Chen, S. H. A., and Lüdecke, D. (2019).
#' Indices of Effect Existence and Significance in the Bayesian Framework.
#' Frontiers in Psychology, 10, 2767. \doi{10.3389/fpsyg.2019.02767}
#'
#' - Rafi Z, Greenland S. Semantic and cognitive tools to aid statistical
#' science: replace confidence and significance by compatibility and surprise.
#' BMC Medical Research Methodology (2020) 20:244.
#'
#' - Schweder T. Confidence is epistemic probability for empirical science.
#' Journal of Statistical Planning and Inference (2018) 195:116–125.
#' \doi{10.1016/j.jspi.2017.09.016}
#'
#' - Schweder T, Hjort NL. Frequentist analogues of priors and posteriors.
#' In Stigum, B. (ed.), Econometrics and the Philosophy of Economics: Theory
#' Data Confrontation in Economics, pp. 285-217. Princeton University Press,
#' Princeton, NJ, 2003
#'
#' - Vos P, Holbert D. Frequentist statistical inference without repeated sampling.
#' Synthese 200, 89 (2022). \doi{10.1007/s11229-022-03560-x}
#'
#' @return A data frame.
#'
#' @examplesIf requireNamespace("bayestestR") && require("see", quietly = TRUE)
#' data(qol_cancer)
#' model <- lm(QoL ~ time + age + education, data = qol_cancer)
#' p_direction(model)
#'
#' result <- p_direction(model)
#' plot(result)
#' @export
p_direction.lm <- function(x,
ci = 0.95,
method = "direct",
null = 0,
...) {
# first, we need CIs
out <- ci(x, ci = ci, ...)
# we now iterate all confidence intervals and create an approximate normal
# distribution that covers the CI-range.
posterior <- as.data.frame(lapply(seq_len(nrow(out)), function(i) {
ci_range <- as.numeric(out[i, c("CI_low", "CI_high")])
.generate_posterior_from_ci(ci, ci_range)
}))
colnames(posterior) <- out$Parameter

out$pd <- as.numeric(bayestestR::p_direction(
posterior,
method = method,
null = null,
...
))

# check we don't have duplicated columns in "posterior" we need this for
# plotting
if (anyDuplicated(colnames(posterior)) > 0 && !is.null(out$Component)) {
comps <- .rename_values(out$Component, "zero_inflated", "zi")
comps <- .rename_values(comps, "conditional", "cond")
colnames(posterior) <- paste0(out$Parameter, "_", comps)
out$Parameter <- paste0(out$Parameter, "_", comps)
}

# we need these for plotting
if (!"Effects" %in% colnames(out)) {
out$Effects <- "fixed"
}
if (!"Component" %in% colnames(out)) {
out$Component <- "conditional"
}

# reorder
out <- out[intersect(c("Parameter", "CI", "CI_low", "CI_high", "pd", "Effects", "Component"), colnames(out))]

attr(out, "data") <- posterior
attr(out, "null") <- null
class(out) <- c("p_direction_lm", "p_direction", "see_p_direction", "data.frame")
out
}


# methods ---------------------------------------------------------------------

#' @export
print.p_direction_lm <- function(x, digits = 2, p_digits = 3, ...) {
null <- attributes(x)$null
caption <- sprintf(
"Probability of Direction (null: %s)",
insight::format_value(null, digits = digits, protect_integer = TRUE)
)

# deal with Effects and Component columns
if ("Effects" %in% colnames(x) && insight::n_unique(x$Effects) == 1) {
x$Effects <- NULL
}
if ("Component" %in% colnames(x) && insight::n_unique(x$Component) == 1) {
x$Component <- NULL
}

x <- insight::format_table(x, digits = digits, p_digits = p_digits)
cat(insight::export_table(x, title = caption, ...))
}


# other classes --------------------------------------------------------------

#' @export
p_direction.glm <- p_direction.lm

#' @export
p_direction.coxph <- p_direction.lm

#' @export
p_direction.svyglm <- p_direction.lm

#' @export
p_direction.glmmTMB <- p_direction.lm

#' @export
p_direction.merMod <- p_direction.lm

#' @export
p_direction.wbm <- p_direction.lm

#' @export
p_direction.lme <- p_direction.lm

#' @export
p_direction.gee <- p_direction.lm

#' @export
p_direction.gls <- p_direction.lm

#' @export
p_direction.feis <- p_direction.lm

#' @export
p_direction.felm <- p_direction.lm

#' @export
p_direction.mixed <- p_direction.lm

#' @export
p_direction.hurdle <- p_direction.lm

#' @export
p_direction.zeroinfl <- p_direction.lm

#' @export
p_direction.rma <- p_direction.lm
49 changes: 37 additions & 12 deletions R/p_significance.R
Original file line number Diff line number Diff line change
Expand Up @@ -12,14 +12,17 @@ bayestestR::p_significance
#' negligible effect in the median's direction, considering a parameter's _full_
#' confidence interval. In other words, it returns the probability of a clear
#' direction of an effect, which is larger than the smalles effect size of
#' interest (e.g., a minimal important difference). In comparison the the
#' [`equivalence_test()`] function, where the *SGPV* (second generation p-value)
#' describes the proportion of the _full_ confidence interval that is _inside_
#' the ROPE, the value returned by `p_significance()` describes the _larger_
#' proportion of the _full_ confidence interval that is _outside_ the ROPE. This
#' makes `p_significance()` comparable to [`bayestestR::p_direction()`],
#' however, while `p_direction()` compares to a point-null by default,
#' `p_significance()` compares to a range-null.
#' interest (e.g., a minimal important difference). Its theoeretical range is
#' from zero to one, but the *ps* is typically larger than 0.5 (to indicate
#' practical significance).
#'
#' In comparison the the [`equivalence_test()`] function, where the *SGPV*
#' (second generation p-value) describes the proportion of the _full_ confidence
#' interval that is _inside_ the ROPE, the value returned by `p_significance()`
#' describes the _larger_ proportion of the _full_ confidence interval that is
#' _outside_ the ROPE. This makes `p_significance()` comparable to
#' [`bayestestR::p_direction()`], however, while `p_direction()` compares to a
#' point-null by default, `p_significance()` compares to a range-null.
#'
#' @param x A statistical model.
#' @inheritParams bayestestR::p_significance
Expand Down Expand Up @@ -118,7 +121,9 @@ bayestestR::p_significance
#' - Vos P, Holbert D. Frequentist statistical inference without repeated sampling.
#' Synthese 200, 89 (2022). \doi{10.1007/s11229-022-03560-x}
#'
#' @return A data frame.
#' @return A data frame with columns for the parameter names, the confidence
#' intervals and the values for practical significance. Higher values indicate
#' more practical significance (upper bound is one).
#'
#' @examplesIf requireNamespace("bayestestR") && packageVersion("bayestestR") > "0.14.0"
#' data(qol_cancer)
Expand All @@ -144,6 +149,23 @@ p_significance.lm <- function(x, threshold = "default", ci = 0.95, verbose = TRU
}))
colnames(posterior) <- out$Parameter

# deal with Effects and Component columns
if ("Effects" %in% colnames(out) && insight::n_unique(out$Effects) == 1) {
out$Effects <- NULL
}
if ("Component" %in% colnames(out) && insight::n_unique(out$Component) == 1) {
out$Component <- NULL
}

# check we don't have duplicated columns in "posterior" we need this for
# plotting
if (anyDuplicated(colnames(posterior)) > 0 && !is.null(out$Component)) {
comps <- .rename_values(out$Component, "zero_inflated", "zi")
comps <- .rename_values(comps, "conditional", "cond")
colnames(posterior) <- paste0(out$Parameter, "_", comps)
out$Parameter <- paste0(out$Parameter, "_", comps)
}

# calculate the ROPE range
if (all(threshold == "default")) {
threshold <- bayestestR::rope_range(x, verbose = verbose)
Expand All @@ -160,6 +182,9 @@ p_significance.lm <- function(x, threshold = "default", ci = 0.95, verbose = TRU
threshold <- 0.1
}

# reorder
out <- out[intersect(c("Parameter", "CI", "CI_low", "CI_high", "ps", "Effects", "Component"), colnames(out))]

attr(out, "data") <- posterior
attr(out, "threshold") <- threshold
class(out) <- c("p_significance_lm", "p_significance", "see_p_significance", "data.frame")
Expand All @@ -170,7 +195,7 @@ p_significance.lm <- function(x, threshold = "default", ci = 0.95, verbose = TRU
# methods ---------------------------------------------------------------------

#' @export
print.p_significance_lm <- function(x, digits = 2, p_digits = 3, ...) {
print.p_significance_lm <- function(x, digits = 2, ...) {
threshold <- attributes(x)$threshold
# make sure it's numeric
if (!is.numeric(threshold)) {
Expand All @@ -184,8 +209,8 @@ print.p_significance_lm <- function(x, digits = 2, p_digits = 3, ...) {
"Practical Significance (threshold: %s)",
toString(insight::format_value(threshold, digits = 2))
)
x$ps <- insight::format_p(x$ps, name = "ps", digits = p_digits)
x <- insight::format_table(x, digits = digits, p_digits = p_digits)
x$ps <- insight::format_pd(x$ps, name = NULL)
x <- insight::format_table(x, digits = digits)
cat(insight::export_table(x, title = caption, ...))
}

Expand Down
4 changes: 3 additions & 1 deletion man/equivalence_test.lm.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

Loading

0 comments on commit 0d60d99

Please sign in to comment.