---
title: "An Introduction to t_TOST"
subtitle: "A new function for TOST with t-tests"
author: "Aaron R. Caldwell"
date: "`r Sys.Date()`"
output:
  rmarkdown::html_vignette:
    toc: True
vignette: >
  %\VignetteIndexEntry{An Introduction to t_TOST}
  %\VignetteEncoding{UTF-8}
  %\VignetteEngine{knitr::rmarkdown}
editor_options:
  chunk_output_type: console
bibliography: references.bib
---

```{r setup, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
library(TOSTER)
library(ggplot2)
library(ggdist)

```

In an effort to make `TOSTER` more informative and easier to use, I created the functions `t_TOST` and `simple_htest`. These function operates very similarly to base R's `t.test` function with a few exceptions. First, `t_TOST` performs 3 t-tests (one two-tailed and two one-tailed tests). Second, `simple_htest` allows you to run equivalence testing or minimal effects testing using a t-test or Wilcoxon-Mann-Whitney tests using the `alternative` argument and the output is the same as `t.test` or `wilcox.test` (in that the object is of the class `htest`). In addition, these functions have a generic method where two vectors can be supplied or a formula can be given (e.g.,`y ~ group`). These functions make it easier to switch between types of t-tests. All three types (two sample, one sample, and paired samples) can be performed/calculated from the same function. Moreover, the summary information and visualizations have been upgraded. This should make the decisions derived from the function more informative and user-friendly. 

These functions are not limited to equivalence tests. Minimal effects testing (MET) is possible. MET is useful for situations where the hypothesis is about a minimal effect and the null hypothesis *is* equivalence.

```{r tostplots,echo=FALSE, message = FALSE, warning = FALSE, fig.show='hold'}

ggplot() +
  geom_vline(aes(xintercept = -.5),
             linetype = "dashed") +
  geom_vline(aes(xintercept = .5),
             linetype = "dashed") +
  geom_text(aes(
    y = 1,
    x = -0.5,
    vjust = -.9,
    hjust = "middle"
  ),
  angle = 90,
  label = 'Lower Bound') +
  geom_text(aes(
    y = 1,
    x = 0.5,
    vjust = 1.5,
    hjust = "middle"
  ),
  angle = 90,
  label = 'Upper Bound') +
  geom_text(aes(
    y = 1,
    x = 0,
    vjust = 1.5,
    hjust = "middle"
  ),
  #alignment = "center",
  label = "H0"
  ) +
  geom_text(aes(
    y = 1,
    x = 1.5,
    vjust = 1.5,
    hjust = "middle"
  ),
  #alignment = "center",
  label = "H1"
  ) +
  geom_text(aes(
    y = 1,
    x = -1.5,
    vjust = 1.5,
    hjust = "middle"
  ),
  #alignment = "center",
  label = "H1"
  ) +
theme_tidybayes() +
  scale_y_continuous(limits = c(0,1.75)) +
  scale_x_continuous(limits = c(-2,2)) +
  labs(x = "", y = "",
       title="Minimal Effect Test",
       caption = "H1 = Alternative Hypothesis \n H0 = Null Hypothesis") +
  theme(
    strip.text = element_text(face = "bold", size = 10),
    axis.text.y = element_blank(),
    axis.ticks.y = element_blank()
  )

ggplot() +
  geom_vline(aes(xintercept = -.5),
             linetype = "dashed") +
  geom_vline(aes(xintercept = .5),
             linetype = "dashed") +
  geom_text(aes(
    y = 1,
    x = -0.5,
    vjust = -.9,
    hjust = "middle"
  ),
  angle = 90,
  label = 'Lower Bound') +
  geom_text(aes(
    y = 1,
    x = 0.5,
    vjust = 1.5,
    hjust = "middle"
  ),
  angle = 90,
  label = 'Upper Bound') +
  geom_text(aes(
    y = 1,
    x = 0,
    vjust = 1.5,
    hjust = "middle"
  ),
  #alignment = "center",
  label = "H1"
  ) +
  geom_text(aes(
    y = 1,
    x = 1.5,
    vjust = 1.5,
    hjust = "middle"
  ),
  #alignment = "center",
  label = "H0"
  ) +
  geom_text(aes(
    y = 1,
    x = -1.5,
    vjust = 1.5,
    hjust = "middle"
  ),
  #alignment = "center",
  label = "H0"
  ) +
theme_tidybayes() +
  scale_y_continuous(limits = c(0,1.75)) +
  scale_x_continuous(limits = c(-2,2)) +
  labs(x = "",
       y = "",
       title="Equivalence Test",
       caption = "H1 = Alternative Hypothesis \n H0 = Null Hypothesis") +
  theme(
    strip.text = element_text(face = "bold", size = 10),
    axis.text.y = element_blank(),
    axis.ticks.y = element_blank()
  )
```

In the general introduction to this package, we detailed how to look at *old* results and how to apply TOST to interpreting those results. However, in many cases, users may have new data that needs to be analyzed. Therefore, `t_TOST` and `simple_htest` can be applied to new data. This vignette will use the `iris` and the `sleep` data.

```{r}
data('sleep')
data('iris')
```

# Independent Groups

For this example, we will use the sleep data. In this data there is a `group` variable and an outcome `extra`.

```{r}
head(sleep)
```

We will assume the data are independent, and that we have equivalence bounds of +/- 0.5 raw units. All we need to do is provide the `formula`, `data`, and  `eqb` arguments for the function to run appropriately. In addition, we can set the `var.equal` argument (to assume equal variance), and the `paired` argument (sets if the data is paired or not). Both are logical indicators that can be set to TRUE or FALSE. The `alpha` is automatically set to 0.05 but this can also be adjusted by the user. The Hedges correction is also automatically calculated, but this can be overridden with the `bias_correction` argument. The `hypothesis` is automatically set to "EQU" for equivalence but if a minimal effect is of interest then "MET" can be supplied. Note: for this example, we will set `smd_ci` to "goulet" since it will reduce the time to produce plots.

```{r}
res1 = t_TOST(formula = extra ~ group,
              data = sleep,
              eqb = .5,
              smd_ci = "goulet")

res1a = t_TOST(x = subset(sleep,group==1)$extra,
               y = subset(sleep,group==2)$extra,
               eqb = .5)
```

We can also using the "simpler" approach with `simple_htest`.

```{r}
# Simple htest

res1b = simple_htest(formula = extra ~ group,
                     data = sleep,
                     mu = .5, # set equivalence bound
                     alternative = "e")


```

Once the function has run, we can print the results with the `print` command. This provides a verbose summary of the results.

```{r}

# t_TOST
print(res1)

# htest

print(res1b)
```


## Plots

Another nice feature is the generic `plot` method that can provide a visual summary of the results (only available for `t_TOST`). All of the plots in this package were inspired by the [concurve](https://cran.r-project.org/package=concurve) R package. There are two types of plots that can be produced. The first, and default, is the consonance density plot (`type = "cd"`).

```{r fig.width=6, fig.height=6}
plot(res1, type = "cd")
```

The shading pattern can be modified with the `ci_shades`.

```{r fig.width=6, fig.height=6, eval=FALSE}
# Set to shade only the 90% and 95% CI areas
plot(res1, type = "cd",
     ci_shades = c(.9,.95))
```

Consonance plots, where all confidence intervals can be simultaneous plotted, can also be produced. The advantage here is multiple confidence interval lines can plotted at once.

```{r fig.width=6, fig.height=6}
plot(res1, type = "c",
     ci_lines =  c(.9,.95))
```

## Descriptions

A description of the results can also be produced with the `describe` or `describe_htest` method and function respectively.

```{r, eval = FALSE}
describe(res1)

describe_htest(res1b)
```

> `r describe(res1)`

> `r describe_htest(res1b)`

# Paired Samples

To perform a paired samples TOST, the process does not change much. We could process the test the same way by providing a formula. All we would need to then is change `paired` to TRUE.

```{r}
res2 = t_TOST(formula = extra ~ group,
              data = sleep,
              paired = TRUE,
              eqb = .5)
res2

res2b = simple_htest(
  formula = extra ~ group,
  data = sleep,
  paired = TRUE,
  mu = .5,
  alternative = "e")
res2b
```

However, we may have two vectors of data that are paired. So we may want to just provide those separately rather than using a data set and setting the formula. This can be demonstrated with the "iris" data.

```{r}
res3 = t_TOST(x = iris$Sepal.Length,
              y = iris$Sepal.Width,
              paired = TRUE,
              eqb = 1)
res3

res3a = simple_htest(
  x = iris$Sepal.Length,
  y = iris$Sepal.Width,
  paired = TRUE,
  mu = 1,
  alternative = "e"
)
res3a
```

We may want to perform a Minimal Effect Test with the `hypothesis` argument set to "MET".

```{r}
res_met = t_TOST(x = iris$Sepal.Length,
              y = iris$Sepal.Width,
               paired = TRUE,
               hypothesis = "MET",
               eqb = 1,
              smd_ci = "goulet")
res_met

res_metb = simple_htest(x = iris$Sepal.Length,
                       y = iris$Sepal.Width,
                       paired = TRUE,
                       mu = 1,
                       alternative = "minimal.effect")
res_metb
```

## Descriptions

A description of the results can also be produced with the `describe` or `describe_htest` method and function respectively.

```{r, eval = FALSE}
describe(res_met)

describe_htest(res_metb)
```

> `r describe(res_met)`

> `r describe_htest(res_metb)`

# One Sample t-test

In other cases we may just have a one sample test. If that is the case all we have to do is supply the `x` argument for the data. For this test we may hypothesis that the mean of Sepal.Length is not more than 5.5 points greater or less than 8.5.

```{r}
res4 = t_TOST(x = iris$Sepal.Length,
              hypothesis = "EQU",
              eqb = c(5.5,8.5),
              smd_ci = "goulet")
res4
```


# Only have the summary statistics? No problem!

In some cases you may only have access to the summary statistics. Therefore, we created a function, `tsum_TOST`, to perform the same tests just based on the summary statistics. This involves providing the function with a number of different arguments.

* `n1 & n2` the sample sizes (only n1 needs to be provided for one sample case)
* `m1 & m2` the sample means
* `sd1 & sd2` the sample standard deviation
* `r12` the correlation between the paired samples; only needed if `paired` is set to TRUE

The results from above can be replicated with the `tsum_TOST`

```{r}
res_tsum = tsum_TOST(
  m1 = mean(iris$Sepal.Length, na.rm=TRUE),
  sd1 = sd(iris$Sepal.Length, na.rm=TRUE),
  n1 = length(na.omit(iris$Sepal.Length)),
  hypothesis = "EQU",
  eqb = c(5.5,8.5)
)

res_tsum
```

```{r fig.width=6, fig.height=6}
plot(res_tsum)
```

```{r}
describe(res_tsum)
```


# Power Analysis for t-test based TOST

We also created `power_t_TOST`  to allow for power calculations for TOST analyses that utilize t-tests. This function uses a more accurate method than the older functions in TOSTER and match the results of the commercially available PASS software. The exact calculations of power are based on Owen’s Q-function or by direct integration of the bivariate non-central t-distribution^[Inspired by @PowerTOST in the `PowerTOST` R package. Please see this package for more options]. Approximate power is implemented via the non-central t-distribution or the ‘shifted’ central t-distribution [@phillips1990, @diletti1991]. The function is limited to power analyses involves one sample, two sample, and paired sample cases. More options are available in the `PowerTOST` R package.

The interface for this function is quite simple and was intended to mimic the base R function `power.t.test`. The user must specify the 2 equivalence bounds, and leave only one of the other options blank (`alpha`, `power`, or `n`). The "true difference" can be set with `delta` and the standard deviation (default is 1) can be set with the `sd` argument. Once everything is set and the function is run, a object of the `power.htest` class will be returned.

As an example, let's say we are looking at an equivalence study where we assume the *true* difference is *at least* 1 unit, the standard deviation is 2.5, and we set the equivalence bounds to 2.5 units as well. If we want to find the sample size adequate to have 95% power at an alpha of 0.025 we enter the following:

```{r}
power_t_TOST(n = NULL,
  delta = 1,
  sd = 2.5,
  eqb = 2.5,
  alpha = .025,
  power = .95,
  type = "two.sample")
```

From the analysis above we would conclude that adequate power is achieved with 74 participants per group and 148 participants in total.

# References