1 Sample Size Determination Under Precision-Based Conditions
The bayes_adcock() function determines the assurance
under precision-based conditions discussed in Adcock (1997). Recall in
the frequentist setting, if \(X_i \sim
N(\theta, \sigma^2)\) for \(i = 1,...,
n\) observations and variance \(\sigma^2\) is known, the precision can be
calculated using \(d =
z_{1-\alpha/2}\frac{\sigma}{\sqrt{n}}\), where \(z_{1-\alpha/2}\) is the critical value for
the \(100(1-\alpha/2)\%\) quartile of
the standard normal distribution. Simple rearranging leads to the
following expression for sample size: \[\begin{equation}
n = z_{1-\alpha/2}^2\frac{\sigma^2}{d^2}\;.
\end{equation}\]
Analysis Stage
Given a random sample with mean \(\bar{x}\), suppose the goal is to estimate population mean \(\theta\).
To formulate this problem in the Bayesian setting, suppose \(x_1, \cdots, x_n\) is a random sample from \(N(\theta, \sigma^2)\) and the sample mean is distributed as \(\bar{x}|\theta, \sigma^2 \sim N(\theta, \sigma^2/n)\). The analysis objective is to observe that the absolute difference between \(\bar{x}\) and \(\theta\) falls within a margin of error no greater than \(d\). Given data \(\bar{x}\) and a pre-specified confidence level \(\alpha\), the assurance can be formally expressed as \[\begin{equation} \delta = P_{\bar{x}}\{\bar{x}: P(|\bar{x} - \theta| \leq d) \geq 1-\alpha\}\;. \end{equation}\] We assign \(\theta \sim N(\theta_0^{(a)}, \sigma^2/n_a)\) as the analysis prior, where \(n_a\) quantifies the amount of prior information we have for \(\theta\). The posterior of \(\theta\) is obtained by taking the product of the prior and likelihood, resulting in \[\begin{equation} N\Bigg(\bar{x} {\left| \theta, \frac{\sigma^2}{n}\right.}\Bigg) \times N\Bigg(\theta {\left | \theta_0^{(a)}, \frac{\sigma^2}{n_a}\right.}\Bigg) = N\Bigg(\theta {\left | \lambda, \frac{\sigma^2}{n_a+n}\right.}\Bigg), \end{equation}\] where \(\lambda = \frac{n\bar{x}+n_a\theta_0^{(a)}}{n_a+n}\). From here we can further evaluate the condition using parameters from the posterior of \(\theta\) to obtain a more explicit version of the analysis stage objective. From \(P(|\bar{x} - \theta| \leq d) = P(\bar{x} - d \leq \theta \leq \bar{x} + d)\), we can standardize all components of the inequality using the posterior distribution of \(\theta\), eventually leading to \[\begin{equation} \label{eq:adcock_region_simplified} \delta = \left\{ \bar{x}: \Phi\left[\frac{\sqrt{n_a + n}}{\sigma} (\bar{x} + d - \lambda)\right] - \Phi\left[\frac{\sqrt{n_a+n}}{\sigma}(\bar{x} - d - \lambda)\right] \geq 1-\alpha \right\}. \end{equation}\]
Design Stage
In the design stage, we need to construct a protocol for sampling
data that will be used to evaluate the analysis objective. This is
achieved by first setting a separate design stage prior on \(\theta\) such that \(\theta \sim N(\theta_0^{(d)},
\sigma^2/n_d)\), where \(n_d\)
quantifies our degree of belief towards the population from which the
sample will be drawn. Given that \(\bar{x}|\theta, \sigma^2 \sim N(\theta,
\sigma^2/n)\), the marginal distribution of \(\bar{x}\) can be computed using
straightforward substitution based on \(\bar{x} = \theta + \epsilon;\) \(\epsilon \sim N(0, \sigma^2/n)\) and \(\theta = \theta_0^{(d)} + \omega;\) \(\omega \sim N(0, \sigma^2/n_d)\).
Substitution \(\theta\) into the
expression for \(\bar{x}\) gives us
\(\bar{x}=\theta_0^{(d)} + (\omega +
\epsilon);\) \((\omega + \epsilon) \sim
N\Big(0, \frac{\sigma^2(n_d+n)}{n_dn}\Big) = N(0, \sigma^2/p)\),
where \(1/p = 1/n_d + 1/n\). The
marginal of \(\bar{x}\) is therefore
\(N(\bar{x}|\theta_0^{(d)},
\sigma^2/p)\), where we will be iteratively drawing samples from
to check if they satisfy the condition outlined in the analysis
stage.
1.1 Example
First, load in the bayesassurance package.
library(bayesassurance)Specify the following inputs:
n: sample size (either vector or scalar).d: precision level in in analysis objective \(|\bar{x} - \theta| \leq d\)mu_beta_a: analysis stage meanmu_beta_d: design stage meann_a: sample size at analysis stage. Quantifies the amount of prior information we have for parameter \(\theta\).n_d: sample size at design stage. Quantifies the amount of prior information we have for where the data is being generated from.sig_sq: known variance \(\sigma^2\).alpha: significance levelmc_iter: number of MC samples evaluated under the analysis objective.
As an example, we assign the following set of arbitrary inputs to
pass into bayes_adcock() and save the output as
out. Notice that we assign a vector of values to sample
size n and choose arbitrary values for the remaining
parameters.
library(ggplot2)
n <- seq(20, 145, 5)
set.seed(20)
out <- bayesassurance::bayes_adcock(n = n, d = 0.20, mu_beta_a = 0.64, mu_beta_d = 0.9,
n_a = 20, n_d = 10, sig_sq = 0.265, alpha = 0.05,
mc_iter = 10000)Within out contains a list of three objects:
assurance_table: table of sample sizes and corresponding assurance values.assur_plot: plot of assurance values. Returned if length ofnis greater than 1.mc_iter: number of MC samples generated.
The first six entries of the resulting power table is shown by
calling head(out$assurance_table).
head(out$assurance_table)| n | Assurance |
|---|---|
| 20 | 0.2378 |
| 25 | 0.3009 |
| 30 | 0.3664 |
| 35 | 0.4376 |
| 40 | 0.5267 |
| 45 | 0.5981 |
The plot of assurance points is produced using the
ggplot2 package. It displays the inputted sample sizes on
the x-axis and the resulting assurance values on the y-axis.
out$assurance_plot
1.2 Overlapping Behaviors Between Frequentist and Bayesian Settings
To demonstrate how the frequentist setting is a special case of the Bayesian framework constructed under precision-based conditions discussed above, we refer back to the fundamental sample size formula, \[\begin{equation} n = z_{1-\alpha/2}^2\frac{\sigma^2}{d^2}\;, \end{equation}\] and perform some simple rearranging to obtain \[\begin{equation} \label{eq:freq_cond} 2\Phi\Big[\frac{\sqrt{n}}{\sigma}d\Big]-1 \geq 1-\alpha\;, \end{equation}\] which we denote as the frequentist condition. Comparing this to the assurance formula we had derived before, \[\begin{equation} \delta = \left\{ \bar{x}: \Phi\left[\frac{\sqrt{n_a + n}}{\sigma} (\bar{x} + d - \lambda)\right] - \Phi\left[\frac{\sqrt{n_a+n}}{\sigma} (\bar{x} - d - \lambda)\right] \geq 1-\alpha \right\}. \end{equation}\] we notice that both expressions have \(1-\alpha\) isolated on the right hand side. Taking \(n_a = 0\) in our assurance formula, simplification proceeds as \[\begin{align*} &\delta = \left\{ \bar{x}: \Phi\left[\frac{\sqrt{n_a + n}}{\sigma} (\bar{x} + d - \lambda)\right] - \Phi\left[\frac{\sqrt{n_a+n}}{\sigma}(\bar{x} - d - \lambda)\right] \geq 1-\alpha \right\} \\ &\xrightarrow[]{n_a = 0} \left\{ \bar{x}: \Phi\left[\frac{\sqrt{n}} {\sigma}d\right] - \Phi \left[-\frac{\sqrt{n}}{\sigma}d\right] \geq 1 - \alpha \right\} = \left\{ \bar{x}: 2\Phi \left[\frac{\sqrt{n}}{\sigma}d\right] - 1 \geq 1 - \alpha \right\}. \end{align*}\] Notice that the resulting term on the left hand side is now equivalent to that of the frequentist condition. This suggests that if we let \(\theta\) take on a weak analysis prior, i.e. \(n_a = 0\), we revert back to the frequentist setting in the analysis stage.
We construct a simple function that returns the left hand side of the frequentist condition expressed above.
freq_cond <- function(n, d, sig_sq){
sigma <- sqrt(sig_sq)
lhs <- 2 * pnorm(d / (sigma / sqrt(n))) - 1
return(lhs)
}Minor modifications can be implemented in the
bayes_adcock() function such that only the left-hand term
is returned (function edits not shown). We will denote this modified
function as adcock_lhs(). Assigning \(n_a = 0\), \(d =
0.15\) and arbitrary values to all other inputs, we obtain
left-hand side values for both frequentist and Bayesian conditions,
which we denote as lhs1 and lhs2 respectively.
Plots are subsequently produced using the ggplot2
package.
library(ggplot2)
n <- seq(5, 300, 5)
n_a <- 1e-8
n_d <- n_a
d <- 0.15
set.seed(1)
sig_sq <- runif(1, 0, 1)
mu_beta_d <- runif(1, 0, 1)
mu_beta_a <- runif(1, 0, 1)
lhs1 <- freq_cond(n, d, sig_sq)
lhs2 <- sapply(n, function(i) adcock_lhs(n = i, d, mu_beta_a, mu_beta_d, n_a,
n_d, sig_sq))
df1 <- as.data.frame(cbind(n, lhs1))
df2 <- as.data.frame(cbind(n, lhs2))
p <- ggplot(df1) + geom_line(data = df1, alpha = 0.5, aes(x = n, y = lhs1,
color="Adcock Cond"), lwd = 1.2)
p1 <- p + geom_point(data = df2, alpha = 0.5, aes(x = n, y = lhs2,
color="Bayesian Sim"),lwd=1.2) + ylab("Probability of Meeting Analysis Objective") +
xlab("Sample Size n") + ggtitle("Comparison Between Adcock Conditiion and Bayesian
Simulation Results")
p1 <- p1 + geom_label(
label="d = 0.15",
x=25,
y=0.98,
color = "black", size = 3
)
We repeat this using various precision values to showcase the
overlapping scenario. Notice that the overall assurance values tend to
decrease as the precision is lowered. 
