-
-
Notifications
You must be signed in to change notification settings - Fork 36
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Revise the way we calculate % in ROPE (SGPV) for frequentist models? #849
Comments
Tagging @vincentarelbundock too, maybe this is of interest for you, too? |
We've discussed this elsewhere I think, and concluded that giving "posterior" probabilities estimates from MLE models was a bad idea? Why can't we just have the p-value from the equivalence test? |
Not sure we discussed it in this context, though. But we may decide on dropping the |
If it's so useful, might as well switch to a framework that can provide with posterior probabilities 😉 |
It's all one song... 😉 |
I'm not following after a quick read but will circle back. …perennial Brenton refrain about confidence distributions and frequentist posterior distributions 😜… |
I don't have a strong view.
I read that book that you recommend and liked it a lot. Thanks for that! Although, this kind of interpretation feels very attractive, I will say that my impressionistic sense is that it seems a bit "exotic" in the stats community. For software defaults and labels, I'd lean toward more traditional/conservative language, or at least use labels which make it clear that the suggested interpretation departs somewhat from canon. |
I think labeling is more a side topic here. ;-) ci <- c(-0.76, 2.53)
# rope <- c(-2, 2)
overlap <- c(-0.76, 2)
100 * diff(overlap) / diff(ci)
#> [1] 83.89058 That is, we consider the CI as a uniform distribution with equal density for all values inside this interval. Usually, you would consider the CI as roughly normal distribution, so the % in ROPE of that interval cannot be calculated as if it was uniformly distributed. Thus, my suggestion to simulate a normal distribution (based on the lower/upper limit) and then calculate the proportion of that probability mass that is inside the ROPE range. This would no longer be in line with the proposed calculation of the SGPV. But it would give a more realistic picture of the "% in ROPE". We could, by default, just report the p-value of the equivalence test and then add columns like SGPV and % in ROPE when requested by the user, if these measures are too avantgarde for some of you to be shown by default :-P |
Haha! again, don't listen to me much here. I don't really use this function and don't have a strong view. |
I think this can work... I have implemented a % in ROPE based on a simulated normal distribution of the CI. The "new" results from % inside ROPE are very close to what we would get from a bootstrapped model or a model where we simulate estimates. We now have the SGPV, calculated in the "classical" way as proposed by Lakens, and the % inside ROPE, which assumes a normally distributed confidence interval and where the proportion of the probability mass of this distribution inside the ROPE is shown. See the close results from simulated model and normal model. library(parameters)
data(qol_cancer)
model <- lm(QoL ~ time + age + education + phq4 + hospital, data = qol_cancer)
equivalence_test(model, metric = "all")
#> # TOST-test for Practical Equivalence
#>
#> ROPE: [-1.99 1.99]
#>
#> Parameter | 90% CI | SGPV | % in ROPE | Equivalence | p
#> -----------------------------------------------------------------------------
#> (Intercept) | [56.21, 69.05] | < .001 | 0% | Rejected | > .999
#> time | [-0.18, 2.52] | 0.802 | 85.68% | Undecided | 0.160
#> age | [-0.40, 0.07] | > .999 | 100% | Accepted | < .001
#> education [mid] | [ 3.12, 9.12] | < .001 | 0% | Rejected | 0.988
#> education [high] | [ 3.92, 11.06] | < .001 | 0% | Rejected | 0.994
#> phq4 | [-5.75, -4.70] | < .001 | 0% | Rejected | > .999
#> hospital | [-1.52, 8.95] | 0.335 | 27.05% | Undecided | 0.743
x <- simulate_model(model)
equivalence_test(x, ci = 0.9)
#> # Test for Practical Equivalence
#>
#> ROPE: [-1.99 1.99]
#>
#> Parameter | H0 | inside ROPE | 90% HDI
#> -------------------------------------------------------
#> (Intercept) | Rejected | 0.00 % | [55.86 69.00]
#> time | Undecided | 88.22 % | [-0.19 2.65]
#> age | Accepted | 100.00 % | [-0.38 0.09]
#> educationmid | Rejected | 0.00 % | [ 3.00 9.08]
#> educationhigh | Rejected | 0.00 % | [ 3.96 11.05]
#> phq4 | Rejected | 0.00 % | [-5.72 -4.74]
#> hospital | Undecided | 25.78 % | [-1.37 9.19]
model <- glm(vs ~ wt + mpg + cyl, data = mtcars, family = "binomial")
equivalence_test(model, metric = "all")
#> # TOST-test for Practical Equivalence
#>
#> ROPE: [-0.18 0.18]
#>
#> Parameter | 90% CI | SGPV | % in ROPE | Equivalence | p
#> ------------------------------------------------------------------------
#> (Intercept) | [-14.74, 30.58] | 0.008 | 0.95% | Undecided | 0.991
#> wt | [ -0.59, 5.11] | 0.064 | 3.79% | Undecided | 0.964
#> mpg | [ -0.48, 0.61] | 0.332 | 42.74% | Undecided | 0.593
#> cyl | [ -5.25, -0.34] | < .001 | 1.68% | Rejected | 0.983
x <- simulate_model(model)
equivalence_test(x, ci = 0.9)
#> # Test for Practical Equivalence
#>
#> ROPE: [-0.18 0.18]
#>
#> Parameter | H0 | inside ROPE | 90% HDI
#> ------------------------------------------------------
#> (Intercept) | Undecided | 0.89 % | [-15.52 30.46]
#> wt | Undecided | 4.22 % | [ -0.33 5.15]
#> mpg | Undecided | 44.22 % | [ -0.48 0.63]
#> cyl | Rejected | 0.00 % | [ -5.46 -0.48]
# for predictions
model <- lm(QoL ~ time + age + education + phq4 + hospital, data = qol_cancer)
ggeffects::ggpredict(model, "education") |>
equivalence_test(metric = "all")
#> # TOST-test for Practical Equivalence
#>
#> ROPE: [-1.99 1.99]
#>
#> education | Contrast | 95% CI | SGPV | % in ROPE | Equivalence | p
#> ---------------------------------------------------------------------------------
#> low-mid | -6.12 | [ -9.69, -2.56] | < .001 | 0.42% | Rejected | 0.988
#> low-high | -7.49 | [-11.74, -3.24] | < .001 | 0% | Rejected | 0.994
#> mid-high | -1.37 | [ -4.59, 1.86] | 0.596 | 61.05% | Undecided | 0.353 |
This idea came up after a chat with Isabella on Twitter.
See following example at https://easystats.github.io/parameters/reference/equivalence_test.lm.html
If we look at the parameter
time
, we see a CI of(-0.76, 2.53)
and 83.52% inside ROPE. The latter is equivalent to the SGPV. Both SGPV and % in ROPE are calculated the same way (see GitHub repo from Daniel: https://github.com/Lakens/TOST_vs_SGPV/tree/master/functions, related paper https://open.lnu.se/index.php/metapsychology/article/view/933).However, although we actually want the overlap of the interval (probability mass) with the ROPE, it's calculated (for simplicity) based on the overlap of the range.
What we are actually doing is:
But what we actually want, is:
So my question is, to compute the % in ROPE for frequentist models, should we switch to simulating the CI and calculate the % in ROPE based on the actual probability mass of that interval, than just the range? Wouldn't this be more precise? Or am I missing something here?
@easystats/core-team
The text was updated successfully, but these errors were encountered: