# Optimal Sequential Procedure for Early Detection
of Multiple Side Effects

Jiayue Wang and Ben Boukai
Department of Mathematical Sciences, IUPUI
Indianapolis, Indiana, 46202  Email<EMAIL_ADDRESS>Email<EMAIL_ADDRESS>
###### Abstract

In this paper, we propose an optimal sequential procedure for the early
detection of potential side effects resulting from the administration of some
treatment (e.g. a vaccine, say). The results presented here extend previous
results obtained in [17] who study the single side effect case to the case of
two (or more) side effects. While the sequential procedure we employ,
simultaneously monitors several of the treatment’s side effects, the
$(\alpha,\beta)$-optimal test we propose does not require any information
about the inter-correlation between these potential side effects. However, in
all of the subsequent analyses, including the derivations of the exact
expressions of the Average Sample Number (ASN), the Power function, and the
properties of the post-test (or post-detection) estimators, we accounted
specifically, for the correlation between the potential side effects.

In the real-life application (such as post-marketing surveillance), the number
of available observations is large enough to justify asymptotic analyses of
the sequential procedure (testing and post-detection estimation) properties.
Accordingly, we also derive the consistency and asymptotic normality of our
post-test estimators; results which enable us to also provide (asymptotic,
post-detection) confidence intervals for the probabilities of various side-
effects. Moreover, to compare two specific side effects, their relative risk
plays an important role. We derive the distribution of the estimated relative
risk in the asymptotic framework to provide appropriate inference. To
illustrate the theoretical results presented, we provide two detailed examples
based on the data of side effects on COVID-19 vaccine collected in Nigeria
(see [16]).

Keywords:  Bivariate optimal sequential test; ASN; Asymptotic normality;
Relative risk; COVID-19 side effects.

AMS Classification: 62L10, 62L12

##  1 Introduction

In a previous paper, [17], we discuss the early direction problem of a single
side effect resulting from some treatment applications (e.g. COVID-19
vaccination). We introduced an optimal sequential procedure for such a
scenario which matches the fixed sample size of the optimal
$(\alpha,\beta)$-UMP test applicable for such a circumstance. In the present
paper we are extending that approach for the early, sequential, detection of
two (or more) side effects. Unlike the case of a single side effect, which is
captured in a sequential random walk of a single binary process, this
generalization, is captured by a sequential random walk over a lattice, of a
bivariate binary response whose components are not independent.

In the literature, one may find numerous studies involving, mostly,
treatment’s efficacy as well as the treatment’s safety measure, which may also
be viewed as situations involving bivariate ’response’ (i.e. ’treatment-
response’ and say, ’toxicity’), however, in drastically various designs. Some
of these studies were carried as two-stage testing procedures which terminate
once certain thresholds on the treatment efficacy or realized toxicity level
have been met. Other were built on a more complex quantitative bivariate
response in group sequential settings, see for example, [9] who described
bivariate normal response in group sequential tests.

A group sequential design for bivariate binary response was introduced by [11]
and by [12] who improved on the calculations via importance sampling. On the
other hand, [10] constructed two-stage design for bivariate binary response
utilizing minimax optimization. Later, [13] modified that design to cope with
the trade-off between safety and efficacy. To further reduce the sample-size
requirement, [14] proposed a curtailed two-stage design for two dependent
binary responses. [15] improved on this stage-wise design with control of the
error rates $\alpha$ and $\beta$. We point out that their basic set-up is
aimed to ’reject’ in the null hypothesis that the tested treatment is
ineffective or unsafe against the alternative hypothesis is that the treatment
is effective and safe.

In this paper, we introduce a new purely sequential design for the early
detection of the multiple potential side effects of a treatment. As a case in
point, we highlight the ’treatment’ of a vaccination campaign (e.g. COVID-19
vaccination campaign). Nowadays, the development to create new vaccines is
being rapidly improving as new technologies are devised and implemented (i.e.
mRNA). Further, with the growing need for a quick deployment of the new
vaccines to the populous, the concern of several potential side effects of the
vaccine becomes a crucial problem.

In modeling such a set up, we first consider the case with two potential side
effects. To that end, consider a vaccination campaign of a population of size
$N$ (we will expand in this later). We assume that the vaccine may potentially
cause for two possible side effects which we label as side effect $X$ and side
effect $Y$, each realized with some unknown probability (proportion)
$\theta_{x}$ and $\theta_{y}$, $(0<\theta_{x},\theta_{y}<1)$, respectively. It
is desired to cease the vaccination campaign, if their proportion of the side
effects, $\theta_{x}$ or $\theta_{y}$ are too large (or unacceptable), namely,
if $\theta_{x}>\theta_{x}^{0}$ or $\theta_{y}>\theta_{y}^{0}$, for some
’desired’ nominal proportions of these two side effects, $\theta_{x}^{0}$ and
$\theta_{y}^{0}$; otherwise, the vaccination campaign should be continued.
This problem can formally be stated as a (sequential) test of the statistical
hypotheses:

$H_{0}:\theta_{x}\leq\theta_{x}^{0}\ \text{{\bf{and}}}\
\theta_{y}\leq\theta_{y}^{0}\ \text{against}\ H_{1}:\theta_{x}>\theta_{x}^{0}\
\text{{\bf{or}} }\ \theta_{y}>\theta_{y}^{0}.$ (1)

For this vaccination campaign, ’patients’ are sequentially receiving the
vaccine (a ’treatment’ of sorts) and subsequently are classified according to
whether or not they have exhibited side effect $X$, side effect $Y$ or both.
Clearly for the sequential testing the hypotheses (1), one would terminate the
vaccination campaign as soon as it has exhibited too many ’side effects’;
otherwise, one would continue uninterruptedly the vaccination process as long
as $H_{0}$ is not rejected.

In Section $2$ below we present the curtailed bivariate sequential test we
propose for the above hypotheses (1). This test is based on a specific
stopping rule which hinges on the bivariate binary process. We demonstrate the
optimality of the proposed test which achieves the desired probabilities
$(\alpha,\beta)$ of the Type I and Type II errors, all while meeting the
maximal sample-size specification of the respective Uniformly Most Powerful
(UMP) test. The derivations and expressions of several important quantities,
such as the Average Sample Number (ASN), and the corresponding power function
of the test are provided. These derivations account for the fact that the two
components of the bivariate binary process, namely $X$ and $Y$ are not
independent and hence, the expressions depend on three parameters,
$\theta_{x}$, $\theta_{y}$ as well as on the correlation between $X$ and $Y$,
namely $\rho\equiv\rho_{x,y}$. Since the stopping rule used for the bivariate
binary process involves some different termination scenarios, which can be
illustrated via a random walk over a lattice, the exact calculations of ASN
and the power function are intricate and tedious. However, exploiting some
implied relationships between $\theta_{x}$, $\theta_{y}$, and $\rho$, we are
able to provide some tight bounds for the ASN by simple expressions.

For the statistical inference of our optimal bivariate sequential test, we
propose, in Section $4$, post-test (or post-detection) estimators of
$\theta_{x}$ and $\theta_{y}$ and analyze their properties. We derive the
exact expressions for the expectation of the post-test estimators, again,
considering the various possible termination scenarios. Additionally, we
further study the properties of our post-test estimators in an asymptotic
framework and providing their asymptotic normality. This asymptotic bivariate
normality of the post-test estimators is exploited further to simplify the
power calculation as well as determination of the ASN. In addition, we discuss
in Section $5$ the asymptotic normality of the estimator for the relative risk
of these two side effects, and derive the appropriate confidence interval for
it. Finally, we close the paper with two detailed examples based on the data
of side effects on COVID-19 vaccine collected from questionnaires in Nigeria
(see [16]).

##  2 The Optimal Bivariate Sequential Test

Consider a vaccination campaign of a population of $N$ individuals who are
being vaccinated sequentially. Each vaccinated individual is being observed
for the expression of two possible side effects, labeled here as $X$ and $Y$.
Having inspected the first $n$ vaccinated individuals out of $N$, $(1\leq
n\leq N)$, the results are summarized in the following $2\times 2$ table

Table 1: Contingency table of $n$ vaccinated people classified by side effect $X$ and side effect $Y$ |  | $Y$ |  |
---|---|---|---|---
|  | No | Yes |
$X$ | No | $n_{00}$ | $n_{01}$ |
| Yes | $n_{10}$ | $n_{11}$ | $n_{1+}$
|  |  | $n_{+1}$ | $n$

That is, $n_{10}$ of the $n$ vaccinated people have exhibited side effect $X$
only, $n_{01}$ of the vaccinated people have exhibited side effect $Y$ only,
whereas $n_{11}$ of the vaccinated people have exhibited both side effects.
Clearly, $n_{00}$ of the $n$ vaccinated individuals have exhibited neither of
the side effects. For a given $n$, the distribution of these counts is the
multinomial distribution,

${\underset{{}^{\sim}}{n}}=(n_{00},n_{10},n_{01},n_{11})^{\prime}\sim{\cal
MN}(n,\underset{{}^{\sim}}{p}),$ (2)

where $\underset{{}^{\sim}}{p}:=(p_{00},p_{10},p_{01},p_{11})^{\prime}$ is the
vector of the corresponding probabilities, $\sum_{i}\sum_{j}p_{ij}=1$ along
with $\sum_{i}\sum_{j}n_{ij}=n$. For each individual, we denote the
corresponding classification indicator by
${\underset{{}^{\sim}}{z}}^{(k)}=(z_{00,k},z_{10,k},z_{01,k},z_{11,k})^{\prime}$,
where $\sum_{i}\sum_{j}z_{ij,k}=1$, for $k=1,2,\dots$. Clearly,

${\underset{{}^{\sim}}{z}}^{(k)}\sim{\cal
MN}(1,\underset{{}^{\sim}}{p}),\qquad k=1,2,\dots.$ (3)

In view of the sequential nature of the vaccination campaign (and hence of the
sampling process), the corresponding sequence of these indicators,
${\underset{{}^{\sim}}{z}}^{(1)},{\underset{{}^{\sim}}{z}}^{(2)},\dots$ become
available, one–at–a–time or in batches, and the final number of individuals
which are being utilized for inference and decisions may depend, in some
prescribed fashion, on the data available to the experimenter, namely,

${\underset{{}^{\sim}}{n}}=\sum_{k=1}^{n}{\underset{{}^{\sim}}{z}}^{(k)},\ \ \
\text{and}\ \ \ n_{ij}=\sum_{k=1}^{n}z_{ij,k}.$

Clearly, for each ${\underset{{}^{\sim}}{z}}^{(k)}$, $k=1,\dots,n$, the
marginal distributions of the indicator for side effect $X$ (where
$X_{i}=z_{11,i}+z_{10,i}$) and of the indicator for side effect $Y$ (where
$Y_{j}=z_{11,j}+z_{01,j}$) are both Bernoulli random variables, so that,

$X\sim{\cal B}ern(1,\theta_{x}),\ \ \ \text{with}\ \
\theta_{x}:=p_{11}+p_{10},$

and

$Y\sim{\cal B}ern(1,\theta_{y}),\ \ \ \text{with}\ \
\theta_{y}:=p_{11}+p_{01}.$

Note from the outset that the indicator for side effect $X$ and the indicator
for side effect $Y$ are not independent random variables (for each $k$). In
fact, it can be easily verified (see for example [7]), that

$Cov(X,Y)=p_{11}-\theta_{x}\theta_{y}.$ (4)

In the case where the two side effects are statistically independent (rarely),
the probabilities
$\underset{{}^{\sim}}{p}:=(p_{00},p_{10},p_{01},p_{11})^{\prime}$ are fully
specified by the two marginal probabilities $\theta_{x}$ and $\theta_{y}$.
However, when the two side effects are not independent (more common), the
correlation between them is

$\rho=\frac{Cov(X,Y)}{\sqrt{Var(X)Var(Y)}}=\frac{p_{11}-\theta_{x}\theta_{y}}{\sqrt{\theta_{x}(1-\theta_{x})\theta_{y}(1-\theta_{y})}}.$

Clearly, this correlation $(-1<\rho<1)$, implies some structural restrictions
on the parameters and the parameter space. Denoting by
$\Omega_{x}=\theta_{x}/(1-\theta_{x})$ and
$\Omega_{y}=\theta_{y}/(1-\theta_{y})$, the odds for exhibiting side effects
$X$ and $Y$, respectively, we have the following three restrictions of $\rho$:

  * •

(i)
$p_{11}=\rho\sqrt{\theta_{x}(1-\theta_{x})\theta_{y}(1-\theta_{y})}+\theta_{x}\theta_{y}\geq
0\ \Rightarrow\ \rho\geq-\sqrt{\Omega_{x}\Omega_{y}}$;

  * •

(ii) $\theta_{x}-p_{11}\geq 0\ \Rightarrow\
\rho\leq\sqrt{\frac{\Omega_{x}}{\Omega_{y}}}$;

  * •

(iii) $\theta_{y}-p_{11}\geq 0\ \Rightarrow\
\rho\leq\sqrt{\frac{\Omega_{y}}{\Omega_{x}}}$.

Therefore, we have the following condition on the correlation $\rho$.

Condition A The values of $\theta_{x}\in(0,1)$, $\theta_{x}\in(0,1)$ and the
correlation $\rho$, between the two side effects $X$ and $Y$ are such

$\left\\{-\sqrt{\Omega_{x}\Omega_{y}}\right\\}\leq\rho\leq\min\left\\{\sqrt{\frac{\Omega_{x}}{\Omega_{y}}},\sqrt{\frac{\Omega_{y}}{\Omega_{x}}}\right\\},$

with $\Omega_{x}={\theta_{x}}/({1-\theta_{x}})$ and
$\Omega_{y}={\theta_{y}}/({1-\theta_{y}})$.

Accordingly, throughout this work, we will restrict attention to the
restricted parameter set as defined by

$\Theta_{R}\equiv\left\\{\theta_{x}\in(0,1),\ \theta_{y}\in(0,1)\text{ and
}|\rho|<1\text{ satisfy {\bf{Condition A}}}\right\\}.$ (5)

Let $S_{n}^{x}=X_{1}+X_{2}+\cdots+X_{n}$, denote the number of people
exhibited side effect $X$, and $S_{n}^{y}=Y_{1}+Y_{2}+\cdots+Y_{n}$, denote
the number of people exhibited side effect $Y$. As was stated earlier, large
enough values of $S_{n}^{x}$ or $S_{n}^{y}$ should lead to the termination of
the vaccinated campaign and to some corrective measures (for the patient’s
treatment). Specifically, in Table 1, the values of $S_{n}^{x}$ and
$S_{n}^{y}$ are displayed by $n_{1+}$ and $n_{+1}$, respectively.

Following [17] (who considered the case of a single side effect), we proceed
by obtaining, for given desired probabilities of Type I error and Type II
error, $(\alpha,\beta)$, the optimal fixed sample size, say $N^{*}$, and a
corresponding critical value $k^{*}$ for the construction of a UMP test for
each of the side effects, separately. For instance, in the case of side effect
$X$, suppose we construct the size $\alpha$ UMP test of
$H_{0}:\theta_{x}\leq\theta_{x}^{0}$ against
$H_{1}:\theta_{x}>\theta_{x}^{0}$, which has a Type II error probability
$\beta$ at some $\theta_{x}=\theta_{x}^{1}>\theta_{x}^{0}$. Standard normal
approximation of the distribution of $S_{n}^{x}$ (see A.1 in Appendix below
for details) leads to the calculated values of the optimal $N^{*}_{x}$ and a
corresponding critical test value $k^{*}_{x}$, for the given
$(\alpha,\beta,\theta_{x}^{0},\theta_{x}^{1})$. Similarly, for side effect
$Y$, one can determine the optimal $N^{*}_{y}$ and $k^{*}_{y}$ which
correspond, to $(\alpha,\beta,\theta_{y}^{0},\theta_{y}^{1})$.

By combining these two separate UMP tests, we consider the construction of an
optimal bivariate sequential testing of (1) which stops the sampling process
as soon as $S^{x}_{n}>k^{*}_{x}$ or $S^{y}_{n}>k^{*}_{y}$. Hence, the sampling
(i.e. vaccination) process is to be terminated at a random _stopping time_
(see [4] for definition),

$M=\min\left\\{M_{x},M_{y}\right\\},$

where

$M_{x}:=\inf\\{n>k^{*}_{x}:\ S^{x}_{n}>k^{*}_{x}\\},\quad
M_{y}:=\inf\\{n>k^{*}_{y}:\ S^{y}_{n}>k^{*}_{y}\\}.$ (6)

The corresponding sequential test of (1) can be written as:

$\operatorname{T_{seq}}:=\begin{cases}\text{if}\ S^{x}_{n}=k_{x}^{*}+1\text{
or }S^{y}_{n}=k_{y}^{*}+1\ &\text{stop and reject $H_{0}$, $M=n$;}\\\ \text{if
}S^{x}_{n}\leq k_{x}^{*}\text{ and }S^{y}_{n}\leq k_{y}^{*}\ \ \
&\text{continue the vaccination.}\end{cases}$ (7)

The following Figure 1 illustrates the random walk over a lattice of two side
effects counting process, which is, the pair $(S_{n}^{x},\,S_{n}^{y})$ jointly
defines a random walk over the integer lattice $\\{1,\dots n\\}^{2}$.

Figure 1: The random walk over a
lattice$O$$S_{n}^{x}$$S_{n}^{y}$$p_{10}$$p_{01}$$p_{11}$$p_{00}$$i$$j$$i+1$$j+1$

However, since in the most realistic situations, the daily supply of vaccines
available to the vaccination center is limited to
$N^{*}:=\min\\{N_{x}^{*},N_{y}^{*}\\}$ units per day (say), the sequential
observation (vaccination) process must be terminated once $N^{*}$ has been
reached. Thus, upon termination, the effective random sample size is

$M^{*}:=\min\left\\{M,N^{*}\right\\}.$

Note that this ’curtailed’ sequential test $\operatorname{T^{*}_{seq}}$ can be
written equivalently in terms of the stopping time $M$ as:

$\operatorname{T^{*}_{seq}}:=\begin{cases}\text{if}\ M\leq N^{*}\ \
&\text{stop and reject $H_{0}$}\\\ \text{if}\ M>N^{*}\ \ &\text{do not reject
$H_{0}$}\end{cases}$ (8)

In Figure 2, we illustrate the individual path for $S_{n}^{x}$ and
$S_{n}^{y}$, upon rejection (black) and upon non-rejection (brown).

Figure 2: The curtailed sequential test
$\operatorname{T^{*}_{seq}}$$O$$S_{n}^{x}$$S_{n}^{y}$$k^{*}_{x}+1$$(k^{*}_{x}+1,k^{*}_{y}+1)$$k^{*}_{y}+1$

To study the properties of our bivariate sequential test of (1), we will
consider the power function of the test $\operatorname{T^{*}_{seq}}$,
evaluated at each
$\underset{{}^{\sim}}{\theta}:=(\theta_{x},\theta_{y})^{\prime}\in\Theta_{R}$,

$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta}):=\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(\operatorname{T^{*}_{seq}}\text{
reject }H_{0}),$

and the Average Sample Number
$ASN^{*}(\underset{{}^{\sim}}{\theta}):=E_{\underset{{}^{\sim}}{\theta}}(M^{*})$.

From (7) and (8), it follows immediately that, $\forall\
\underset{{}^{\sim}}{\theta}\in\Theta_{R}$,

$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})=\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(M\leq
N^{*})=1-\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(M>N^{*})\equiv
1-\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(S^{x}_{N^{*}}\leq
k_{x}^{*}\text{ and }S^{y}_{N^{*}}\leq k_{y}^{*}\right).$ (9)

######  Theorem 1.

Let $\operatorname{T^{*}_{seq}}$ be the bivariate curtailed sequential test of
the hypotheses in (1) as given in (8) above and let
$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$ be its power
function. For given $(\alpha,\beta)$, we let $\tilde{\alpha}=\alpha/2$ and
$\tilde{\beta}=\beta$ be the probabilities of the Type I and Type II errors
with corresponding $(N_{x}^{*},k_{x}^{*})$ and $(N^{*}_{y},k^{*}_{y})$ for
each side effect marginal UMP test with
${\underset{{}^{\sim}}{\theta_{0}}}:=(\theta_{x}^{0},\theta_{y}^{0})$ and
$\underset{{}^{\sim}}{\theta_{1}}:=(\theta_{x}^{1},\theta_{y}^{1})$. Then we
have that $\operatorname{T^{*}_{seq}}$ is optimal in the sense that with
$N^{*}=\min\\{N_{x}^{*},N_{y}^{*}\\}$,

$\Pr\left(\text{Type~{}I error of
}\operatorname{T^{*}_{seq}}\right)\leq\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta_{0}})\leq\alpha$

and

$\Pr\left(\text{Type~{}II error of }\operatorname{T^{*}_{seq}}\right)\leq
1-\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta_{1}})\leq\beta.$

Proof In Lemma 2 we provide that the power function
$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$ is
monotonically increasing with respect to $\theta_{x}$ and $\theta_{y}$. Hence
we have that

$\Pr\left(\text{Type~{}I error of
}\operatorname{T^{*}_{seq}}\right)\leq\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta_{0}}),$

and

$\Pr\left(\text{Type~{}II error of }\operatorname{T^{*}_{seq}}\right)\leq
1-\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta_{1}}).$

Suppose $N_{x}^{*}\leq N_{y}^{*}$, we have
$N^{*}=\min\\{N_{x}^{*},N_{y}^{*}\\}\equiv N_{x}^{*}$. Note that since
$N_{x}^{*}\leq N_{y}^{*}$, we also have
$\operatorname{P}_{\theta_{y}}(S_{N^{*}_{y}}^{y}>a)=\operatorname{P}_{\theta_{y}}(\sum_{i=1}^{N^{*}_{y}}Y_{i}>a)\geq\operatorname{P}_{\theta_{y}}(\sum_{i=1}^{N^{*}_{x}}Y_{i}>a)=\operatorname{P}_{\theta_{y}}(S_{N^{*}_{x}}^{y}>a)$,
$\forall\ a\in{\mathbb{R}}$. Accordingly, by (9), it follows that

$\displaystyle\operatorname{P}_{\underset{{}^{\sim}}{\theta_{0}}}\left(M>N^{*}\right)=$
$\displaystyle\operatorname{P}_{\underset{{}^{\sim}}{\theta_{0}}}\left(M_{x}>N^{*}\text{
and }M_{y}>N^{*}\right)$ $\displaystyle=$
$\displaystyle\operatorname{P}_{\underset{{}^{\sim}}{\theta_{0}}}\left(S_{N^{*}}^{x}\leq
k_{x}^{*}\text{ and }S_{N^{*}}^{y}\leq k_{y}^{*}\right)$ $\displaystyle\geq$
$\displaystyle\operatorname{P}_{\theta_{x}^{0}}\left(S_{N^{*}}^{x}\leq
k_{x}^{*}\right)+\operatorname{P}_{\theta_{y}^{0}}\left(S_{N^{*}}^{y}\leq
k_{y}^{*}\right)-1$ $\displaystyle\geq$
$\displaystyle\left(1-\frac{\alpha}{2}\right)+\left(1-\frac{\alpha}{2}\right)-1$
$\displaystyle\geq$ $\displaystyle 1-\alpha.$

Hence, we conclude that

$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta_{0}})=1-\operatorname{P}_{\underset{{}^{\sim}}{\theta_{0}}}\left(M>N^{*}\right)\leq\alpha.$

On the other hand,

$\displaystyle
1-\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta_{1}})=$
$\displaystyle
1-\operatorname{P}_{\underset{{}^{\sim}}{\theta_{1}}}\left(M\leq N^{*}\right)$
$\displaystyle=$
$\displaystyle\operatorname{P}_{\underset{{}^{\sim}}{\theta_{1}}}\left(S_{N^{*}}^{x}\leq
k_{x}^{*}\text{ and }S_{N^{*}}^{y}\leq k_{y}^{*}\right)$ $\displaystyle\leq$
$\displaystyle\operatorname{P}_{\theta_{x}^{1}}\left(S_{N^{*}}^{x}\leq
k_{x}^{*}\right)$ $\displaystyle\leq$ $\displaystyle\beta.$

The proof is similar in the case of $N_{x}^{*}\geq N_{y}^{*}$. This completes
the proof of the Theorem. $\Box$

In Figure 3 below we illustrate the results of Theorem 1 concerning the power
function $\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$ with
respect to the value of $\underset{{}^{\sim}}{\theta}\in\Theta_{R}$.

Figure 3: The contour plot of the power function with respect to $\theta_{x},\
\theta_{y}\in(0,0.3)$ and $\rho=0.1$ for fixed $\alpha=0.05$, $\beta=0.1$,
$\theta_{x}^{0}=0.05$, $\theta_{x}^{1}=0.1$ and $\theta_{y}^{0}=0.1$,
$\theta_{y}^{1}=0.2$ with $N^{*}=121$, $k^{*}_{x}=19$ and $k^{*}_{y}=18$. The
corresponding
$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta_{0}})=0.0321$
and
$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta_{1}})=0.9065$.

##  3 On the ASN

In order to study the efficiency of our proposed optimal bivariate sequential
test in (8), we derive its Average Sample Number (ASN),
$ASN^{*}(\underset{{}^{\sim}}{\theta}):=E_{\underset{{}^{\sim}}{\theta}}(M^{*})$.
Clearly, since $M^{*}=\min\\{M,N^{*}\\}$, we have

$ASN^{*}(\underset{{}^{\sim}}{\theta})=E_{\underset{{}^{\sim}}{\theta}}\left(\min\\{M,N^{*}\\}\right)=N^{*}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M>N^{*}\right)+E_{\underset{{}^{\sim}}{\theta}}\left(M\mathbbm{1}\left[M\leq
N^{*}\right]\right),$

where $M=\min\\{M_{x},M_{y}\\}$ and where $\mathbbm{1}[{\cal A}]$ is the
indicator function of the ’event’ ${\cal A}$. By Lemma 3, we derive that

$\displaystyle ASN^{*}(\underset{{}^{\sim}}{\theta})=$ $\displaystyle
N^{*}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M>N^{*}\right)+\sum_{m=1}^{N^{*}}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right)-N^{*}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
N^{*}+1\right)$ $\displaystyle=$
$\displaystyle\sum_{m=1}^{N^{*}}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right):=I^{*}_{1}.$ (10)

By direct calculation for the term, $I^{*}_{1}$, we have that with
$\overset{{}_{\sim}}{k}:=\min\\{k^{*}_{x},k^{*}_{y}\\}$

$I^{*}_{1}:=N^{*}-\sum_{m=\overset{{}_{\sim}}{k}+1}^{N^{*}-1}(N^{*}-m)\cdot\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M=m\right).$
(11)

In Remark $3$, we provide the exact calculations of the probabilities of
$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(M=m)$ for
$m=\overset{{}_{\sim}}{k}+1,\dots,N^{*}$, under $\Theta_{R}$ in (5), which
once obtained, provide the explicit expression for the evaluation of
$ASN^{*}(\underset{{}^{\sim}}{\theta})$ in (3). These calculations in (11)
utilize the (’negative’ version of) the multinomial distribution of ${\cal
MN}(n,\underset{{}^{\sim}}{p})$ in (2). However, since the exact calculations
of the ASN in (11) are intricate, we introduce below some tight bounds for its
numerical evaluation.

To simplify the expressions, we denote, following (3),

$U_{1}:=E_{\theta_{x}}(M^{*}_{x})=\sum_{m=1}^{N^{*}}\operatorname{P}_{\theta_{x}}\left(M_{x}\geq
m\right),\qquad
U_{2}:=E_{\theta_{y}}(M^{*}_{y})=\sum_{m=1}^{N^{*}}\operatorname{P}_{\theta_{y}}\left(M_{y}\geq
m\right),$

and

$L_{1}:=\sum_{m=1}^{N^{*}}\operatorname{P}_{\theta_{x}}\left(M_{x}\geq
m\right)\operatorname{P}_{\theta_{y}}\left(M_{y}\geq m\right).$

Note that it always holds that

$L_{1}\leq\min\left\\{U_{1},U_{2}\right\\}.$

Also, since $\forall$ $m\in\mathbbm{N}$,

$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right)=\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\min\left\\{M_{x},M_{y}\right\\}\geq
m\right)\leq\min\left\\{\operatorname{P}_{\theta_{x}}\left(M_{x}\geq
m\right),\operatorname{P}_{\theta_{y}}\left(M_{y}\geq m\right)\right\\},$ (12)

we obtain that

$\sum_{m=1}^{N^{*}}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right)\leq\min\left\\{U_{1},U_{2}\right\\}.$

Further, if $M_{x}$ and $M_{y}$ are positively associated random variables,
then $\forall$ $m\in\mathbbm{N}$,

$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right)=\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{x}\geq
m\text{ and }M_{y}\geq
m\right)\geq\operatorname{P}_{\theta_{x}}\left(M_{x}\geq
m\right)\operatorname{P}_{\theta_{y}}\left(M_{y}\geq m\right),$

and therefore,

$\sum_{m=1}^{N^{*}}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right)\geq L_{1}.$ (13)

On the other hand, if $M_{x}$ and $M_{y}$ are negatively associated random
variables, then $\forall$ $m\in\mathbbm{N}$,

$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right)=\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{x}\geq
m\text{ and }M_{y}\geq
m\right)\leq\operatorname{P}_{\theta_{x}}\left(M_{x}\geq
m\right)\operatorname{P}_{\theta_{y}}\left(M_{y}\geq m\right),$

we therefore will have

$\sum_{m=1}^{N^{*}}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right)\leq L_{1}.$ (14)

Accordingly, by (12) and (13), we have the following bounds for the
$ASN^{*}(\underset{{}^{\sim}}{\theta})$

$L_{1}\leq
ASN^{*}(\underset{{}^{\sim}}{\theta})\leq\min\left\\{U_{1},U_{2}\right\\}.$

whenever $M_{x}$ and $M_{y}$ are positively associated random variables (see
Lemma 4). Otherwise, by (14), we have that

$0<ASN^{*}(\underset{{}^{\sim}}{\theta})\leq L_{1}.$

Clearly, if $M_{x}$ and $M_{y}$ are uncorrelated, we will have

$ASN^{*}(\underset{{}^{\sim}}{\theta})=L_{1}.$

Now, since $U_{1}\equiv ASN^{*}(\theta_{x})$ and $U_{2}\equiv
ASN^{*}(\theta_{y})$, we may obtain, by utilizing Theorem $2$ of [17] (see
expression $(20)$ there for details), that

$U_{1}=N^{*}\mathcal{I}_{1-\theta_{x}}\left(N^{*}-k^{*}_{x},k^{*}_{x}+1\right)+\frac{k^{*}_{x}+1}{\theta_{x}}\mathcal{I}_{\theta_{x}}\left(k^{*}_{x}+2,N^{*}-k^{*}_{x}\right),$

where $\mathcal{I}_{\cdot}(\cdot,\cdot)$ denotes the regularized incomplete
beta function111For any $\xi>0,a>0$ and $b>0$,
$\mathcal{I}_{\xi}(a,b):=\int_{0}^{\xi}f(u|a,b)du\equiv\int_{0}^{\xi}\frac{u^{a-1}(1-u)^{b-1}}{b(a,b)}du$,
with $f(u|a,b)$ being the $pdf$ of a ${\cal{B}}eta(a,b)$ random variable. For
a detailed discussion of this function and its relation to binomial
probabilities, see [1], [2]., and similarly,

$U_{2}=N^{*}\mathcal{I}_{1-\theta_{y}}\left(N^{*}-k^{*}_{y},k^{*}_{y}+1\right)+\frac{k^{*}_{y}+1}{\theta_{y}}\mathcal{I}_{\theta_{y}}\left(k^{*}_{y}+2,N^{*}-k^{*}_{y}\right).$

For the direct calculation of the expression $L_{1}$, we must consider two
possible situations:

  * •

(i) $k^{*}_{x}\geq k^{*}_{y}$:

$L_{1}=\mu_{N^{*},x}-\sum_{i=k^{*}_{y}+1}^{k^{*}_{x}}\mathcal{I}_{\theta_{y}}(k^{*}_{y}+1,i-k^{*}_{y})-\sum_{i=k^{*}_{x}+1}^{N^{*}-1}\mathcal{I}_{\theta_{y}}(k^{*}_{y}+1,i-k^{*}_{y})\mathcal{I}_{1-\theta_{x}}(i-k^{*}_{x},k^{*}_{x}+1);$

  * •

(ii) $k^{*}_{x}\leq k^{*}_{y}$:

$L_{1}=\mu_{N^{*},y}-\sum_{i=k^{*}_{x}+1}^{k^{*}_{y}}\mathcal{I}_{\theta_{x}}(k^{*}_{x}+1,i-k^{*}_{x})-\sum_{i=k^{*}_{y}+1}^{N^{*}-1}\mathcal{I}_{\theta_{x}}(k^{*}_{x}+1,i-k^{*}_{x})\mathcal{I}_{1-\theta_{y}}(i-k^{*}_{y},k^{*}_{y}+1).$

Table 2: The values of $ASN^{*}(\underset{{}^{\sim}}{\theta})$ with respect to $\theta_{x}$ and $\theta_{y}$ for fixed $\alpha=0.05$, $\beta=0.1$, $\theta_{x}^{0}=0.05$, $\theta_{x}^{1}=0.1$ and $\theta_{y}^{0}=0.1$, $\theta_{y}^{1}=0.2$ with $N^{*}=121$, $k^{*}_{x}=19$ and $k^{*}_{y}=18$. | $ASN^{*}(\underset{{}^{\sim}}{\theta})$
---|---
| $\left(\theta_{x}^{0},\theta_{y}^{0}\right)$ | $\left(\theta_{x}^{1},\theta_{y}^{0}\right)$ | $\left(\theta_{x}^{0},\theta_{y}^{1}\right)$ | $\left(\theta_{x}^{1},\theta_{y}^{1}\right)$ | $(\theta_{x}^{0},0.25)$ | $(0.25,\theta_{y}^{0})$ | $(0.25,0.25)$
$\rho=0.1$ | upper | 120.6654 | 120.6654 | 93.8602 | 93.8602 | 75.9630 | 79.9251 | 75.9630
| exact | 120.6653 | 120.5080 | 93.8602 | 93.8397 | 75.9630 | 79.9165 | 69.7126
| lower | 120.6653 | 120.5052 | 93.8602 | 93.8282 | 75.9630 | 79.9095 | 69.2791
$\rho=-0.1$ | upper | NA | 120.5052 | 93.8602 | 93.8282 | 75.9630 | 79.9095 | 69.2791
| exact | NA | 120.5035 | 93.8602 | 93.8140 | 75.9630 | 79.8995 | 68.8663
Note: since $\left(\theta_{x}^{0},\theta_{y}^{0}\right)\notin\Theta_{R}$ when
$\rho=-0.1$, ’NA’ is presented.

In Table 2 above we illustrate the calculated bounds for the
$ASN^{*}(\underset{{}^{\sim}}{\theta})$ in comparison to its exact calculated
value.

After deriving the expression of $ASN^{*}(\underset{{}^{\sim}}{\theta})\equiv
E_{\underset{{}^{\sim}}{\theta}}(M^{*})$, we now consider the second moment
$E_{\underset{{}^{\sim}}{\theta}}\left({M^{*}}^{2}\right)$ in order to derive
$Var_{\underset{{}^{\sim}}{\theta}}\left(M^{*}\right)$.

According Lemma 5, we have that

$\displaystyle E_{\underset{{}^{\sim}}{\theta}}\left({M^{*}}^{2}\right)=$
$\displaystyle{N^{*}}^{2}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M>N^{*}\right)+E_{\underset{{}^{\sim}}{\theta}}\left(M^{2}\mathbbm{1}\left[M\leq
N^{*}\right]\right)$ $\displaystyle=$
$\displaystyle{N^{*}}^{2}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M>N^{*}\right)+2\sum_{m=1}^{N^{*}}\left(m-\frac{1}{2}\right)\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right)-{N^{*}}^{2}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
N^{*}+1\right)$ $\displaystyle=$ $\displaystyle
2\sum_{m=1}^{N^{*}}\left(m-\frac{1}{2}\right)\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\geq
m\right):=I^{*}_{2}.$ (15)

For the exact calculation of the term, $I^{*}_{2}$ above, by direct
derivation, we have that

$I^{*}_{2}:={N^{*}}^{2}-\sum_{m=\overset{{}_{\sim}}{k}+1}^{N^{*}-1}\left({N^{*}}^{2}-m^{2}\right)\cdot\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M=m\right).$
(16)

Therefore, we may obtain the exact value of
$Var_{\underset{{}^{\sim}}{\theta}}(M^{*})$ by combining the value of
$E_{\underset{{}^{\sim}}{\theta}}(M^{*})$ in (3) and the value of
$E_{\underset{{}^{\sim}}{\theta}}({M^{*}}^{2})$ in (3), namely,

$Var_{\underset{{}^{\sim}}{\theta}}\left(M^{*}\right)=E_{\underset{{}^{\sim}}{\theta}}\left({M^{*}}^{2}\right)-\left[E_{\underset{{}^{\sim}}{\theta}}\left(M^{*}\right)\right]^{2}=I_{2}^{*}-\left(I_{1}^{*}\right)^{2}.$

Furthermore, to show the relative dispersion of the terminal sample size,
$M^{*}$, we calculate the value of the coefficient of variation (CV) as

$CV_{\underset{{}^{\sim}}{\theta}}\left(M^{*}\right):=\frac{\sqrt{Var_{\underset{{}^{\sim}}{\theta}}\left(M^{*}\right)}}{E_{\underset{{}^{\sim}}{\theta}}\left(M^{*}\right)}.$

As it appears from Table 3, the proposed optimal bivariate sequential test
results are with relatively low CV with respect to the values of $\theta_{x}$
and $\theta_{y}$.

Table 3: The values of the variance and coefficient of variation (CV) of $M^{*}$ with respect to $\theta_{x}$ and $\theta_{y}$ for fixed $\alpha=0.05$, $\beta=0.1$, $\rho=0.1$, $\theta_{x}^{0}=0.05$, $\theta_{x}^{1}=0.1$ and $\theta_{y}^{0}=0.1$, $\theta_{y}^{1}=0.2$ with $N^{*}=121$, $k^{*}_{x}=19$ and $k^{*}_{y}=18$. $\theta_{x}$ $\theta_{y}$ Variance (CV) | 0.02 | $0.1^{*}$ | $0.2^{*}$ | 0.25 | 0.4
---|---|---|---|---|---
0.02 | 3.7e-10 | 6.1436 | 295.2041 | 224.175 | 71.2500
| (1.6e-07) | (0.0205) | (0.1831) | (0.1971) | (0.1777)
$0.05^{*}$ | 0.0002 | 6.1438 | 295.2041 | 224.1750 | 71.2500
| (0.0001) | (0.0205) | (0.1831) | (0.1971) | (0.1777)
$0.1^{*}$ | 2.7969 | 8.7733 | 294.6476 | 224.0743 | 71.2500
| (0.0138) | (0.0246) | (0.1829) | (0.1971) | (0.1777)
0.25 | 232.7980 | 232.4461 | 176.0087 | 139.2098 | 68.8701
| (0.1909) | (0.1908) | (0.1738) | (0.1692) | (0.1751)
0.4 | 75.0000 | 75.0000 | 74.0820 | 69.6496 | 43.5088
| (0.1732) | (0.1732) | (0.1723) | (0.1680) | (0.1496)
Note: ∗ denotes values under the hypotheses $H_{0}$ and $H_{1}$.

##  4 Estimation

Once that the optimal sequential test of the hypotheses (1) is implemented,
one might be interested in the estimation, upon termination, of the unknown
parameters $\theta_{x}$ and $\theta_{y}$. In this section, we derive the post-
test (or post-detection) estimators $\hat{\theta}_{x}$ and $\hat{\theta}_{y}$
of the two parameters $\theta_{x}$ and $\theta_{y}$ and study their properties
(in the finite sample sense as well as in a suitable asymptotic framework).

We begin with the post-test estimator of $\theta_{x}$, as the derivations of
the estimator of $\theta_{y}$ are similar. Clearly, with the curtailed
’stopping time’, $M^{*}:=\min\\{N^{*},M\\}$ of the sequential test, we
consider the ’sample proportion’ of those who have exhibited side effect $X$,
namely,

$\hat{\theta}_{x}=\frac{S_{M^{*}}^{x}}{M^{*}}\equiv\frac{S_{N^{*}}^{x}}{N^{*}}\mathbbm{1}\left[M>N^{*}\right]+\frac{S_{M}^{x}}{M}\mathbbm{1}\left[M\leq
N^{*}\right].$ (17)

The expectation of $\hat{\theta}_{x}$ is

$E_{\underset{{}^{\sim}}{\theta}}\left(\hat{\theta}_{x}\right)=E_{\underset{{}^{\sim}}{\theta}}\left(\frac{S_{N^{*}}^{x}}{N^{*}}\mathbbm{1}\left[M>{N^{*}}\right]\right)+E_{\underset{{}^{\sim}}{\theta}}\left(\frac{S_{M}^{x}}{M}\mathbbm{1}\left[M\leq
N^{*}\right]\right).$ (18)

The first term in (18) can be directly calculated by utilizing the multinomial
distribution of ${\cal MN}(N^{*},\underset{{}^{\sim}}{p})$ in (2) as

$E_{\underset{{}^{\sim}}{\theta}}\left(\frac{S_{N^{*}}^{x}}{N^{*}}\mathbbm{1}\left[M>N^{*}\right]\right)=\sum_{i=0}^{k^{*}_{x}-z}\sum_{j=0}^{k^{*}_{y}-z}\sum_{z=0}^{\overset{{}_{\sim}}{k}}\frac{z+i}{N^{*}}\frac{N^{*}!}{z!i!j!(N^{*}-z-i-j)!}p_{11}^{z}p_{10}^{i}p_{01}^{j}p_{00}^{N^{*}-z-i-j},$
(19)

The second term in (18) is obtained by utilizing the (’negative’ version of)
the multinomial distribution of ${\cal MN}(n,\underset{{}^{\sim}}{p})$ in (2)
as

$E_{\underset{{}^{\sim}}{\theta}}\left(\frac{S_{M}^{x}}{M}\mathbbm{1}\left[M\leq
N^{*}\right]\right)\\\
=\sum_{m=k^{*}_{x}+1}^{N^{*}}\sum_{i=1}^{\min\\{k^{*}_{x}+1,k^{*}_{y}\\}}\sum_{j=0}^{\min\\{k^{*}_{y}-i,m-k^{*}_{x}-1\\}}\frac{k^{*}_{x}+1}{m}A^{(1)}_{i,j}+\sum_{m=k^{*}_{x}+1}^{N^{*}}\sum_{i=0}^{\overset{{}_{\sim}}{k}}\sum_{j=0}^{\min\\{k^{*}_{y}-i,m-k^{*}_{x}-1\\}}\frac{k^{*}_{x}+1}{m}A^{(2)}_{i,j}\\\
+\sum_{m=k^{*}_{y}+1}^{N^{*}}\sum_{i=1}^{\min\\{k^{*}_{x},k^{*}_{y}+1\\}}\sum_{j=0}^{\min\\{k^{*}_{x}-i,m-k^{*}_{y}-1\\}}\frac{i+j}{m}A^{(3)}_{i,j}+\sum_{m=k^{*}_{y}+1}^{N^{*}}\sum_{i=0}^{\overset{{}_{\sim}}{k}}\sum_{j=0}^{\min\\{k^{*}_{x}-i,m-k^{*}_{y}-1\\}}\frac{i+j}{m}A^{(4)}_{i,j}\\\
+\sum_{i=1}^{\overset{{}_{\sim}}{k}+1}\sum_{m=k^{*}_{x}+k^{*}_{y}+2-i}^{N^{*}}\frac{k^{*}_{x}+1}{m}A^{(5)}_{i},$
(20)

where the detailed expressions of $A^{(1)}_{i,j},\ A^{(2)}_{i,j},\
A^{(3)}_{i,j},\ A^{(4)}_{i,j},\ A^{(5)}_{i}$ are given in Remark 3.

Therefore, combining (19) and (20) together, we obtain the exact value of
expectation of the post-test estimator $\hat{\theta}_{x}$ in (18). In Section
$5$, we discuss the performance of our post-test estimators by way of the
maximal relative absolute bias and its tendency between the estimated and the
actual value.

##  5 Asymptotic Properties

We study the asymptotic behaviors of our optimal bivariate sequential test in
a ’local asymptotic’ sense, as $\theta_{x}^{1}\to\theta_{x}^{0}$ and
$\theta_{y}^{1}\to\theta_{y}^{0}$. Specifically, as
$\theta_{x}^{1}\equiv\theta_{x}^{0}(1+\delta)$ and
$\theta_{y}^{1}\equiv\theta_{y}^{0}(1+\delta)$ for some $\delta>0$ and small,
with $\delta\rightarrow 0$222Alternatively, one may take
$\theta_{x}^{1}=\theta_{x}^{0}+\delta$ and
$\theta_{y}^{1}=\theta_{y}^{0}+\delta$, in the subsequent derivations which
lead to the same asymptotic results described here.. Under this
parameterization, we denote by $N_{x}^{\delta}\equiv
N^{*}_{x}(\tilde{\alpha},\beta,\theta_{x}^{0},\delta)$ and
$k_{x}^{\delta}\equiv k_{x}^{*}(\tilde{\alpha},\beta,\theta_{x}^{0},\delta)$,
the optimal $N$ and $k$ of the marginal test on the side effect $X$ and and
similarly, by $N_{y}^{\delta}$ and $k_{y}^{\delta}$ for the marginal test on
the side effect $Y$.

By Lemma $1$ of [17], which discussed the case of the single side effect, we
have as $\delta\to 0$

$N_{x}^{\delta}\to\infty,\quad k_{x}^{\delta}\to\infty,\quad
N_{y}^{\delta}\to\infty,\quad k_{y}^{\delta}\to\infty.$

Similarly, we denote by
$\overset{{}_{\sim}}{k}_{\delta}:=\min\\{k_{x}^{\delta},k_{y}^{\delta}\\}$.
Therefore, we immediately have as $\delta\to 0$,

$N^{*}_{\delta}=\min\left\\{N_{x}^{\delta},N_{y}^{\delta}\right\\}\to\infty\quad\text{and}\quad
M^{*}_{\delta}\geq\overset{{}_{\sim}}{k}_{\delta}+1\to\infty.$

To approximate the power function
$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$ in (9), we
note at first that the power function
$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$ can be
written in an equivalent form as

$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})\equiv
1-\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{\delta}>N_{\delta}^{*}\right)\equiv
1-\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(S_{N^{*}_{\delta}}^{x}\leq
k_{x}^{\delta}\ \text{and}\ \ S_{N^{*}_{\delta}}^{y}\leq
k_{y}^{\delta}\right).$

Further, since $N^{*}_{\delta}\to\infty$ as $\delta\to 0$, we may utilize the
multivariate version of the CLT, along with lemma 2. It is straightforward to
verify that as $\delta\to 0$,

$\begin{pmatrix}[1.5]S_{N^{*}_{\delta}}^{x}\\\
S_{N^{*}_{\delta}}^{y}\end{pmatrix}\sim\mathcal{N}_{2}\left(\bm{\mu_{1}},\bm{V_{1}}\right),$
(21)

where $\mathcal{N}_{2}\left(\bm{\mu_{1}},\bm{V_{1}}\right)$ denotes that
bivariate normal distribution whose mean and variance-covariance matrix are,

$\bm{\mu_{1}}=\begin{pmatrix}[1.5]N^{*}_{\delta}\theta_{x}\\\
N^{*}_{\delta}\theta_{y}\end{pmatrix},\quad\bm{V_{1}}=\begin{pmatrix}[1.5]N^{*}_{\delta}\theta_{x}(1-\theta_{x})&N^{*}_{\delta}\rho\sqrt{\theta_{x}(1-\theta_{x})\theta_{y}(1-\theta_{y})}\\\
N^{*}_{\delta}\rho\sqrt{\theta_{x}(1-\theta_{x})\theta_{y}(1-\theta_{y})}&N^{*}_{\delta}\theta_{y}(1-\theta_{y})\end{pmatrix}.$

Accordingly, the power function
$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$ can be
approximated (when incorporated the standard continuity correction), as
$\delta\to 0$, by

$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})=1-\int_{-\infty}^{k_{x}^{\delta}+0.5}\int_{-\infty}^{k_{y}^{\delta}+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{1}},\bm{V_{1}}\right)dwdu,$
(22)

where $\phi_{2}(u,w\mid\bm{\mu_{1}},\bm{V_{1}})$ denotes the pdf of the
bivariate normal distribution
$\mathcal{N}_{2}\left(\bm{\mu_{1}},\bm{V_{1}}\right)$, above.

The next lemma is a restatement of Theorem $3$ (i) and Theorem $4$ of [5],
which is critical for the derivations of the asymptotic approximation to the
$pmf$ of the stopping time $M_{\delta}$, namely of
$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(M_{\delta}=m)$.

As defined above, let
$\left\\{\left(X_{i},Y_{i}\right)\right\\}_{i=1}^{\infty}$ be i.i.d. two-
dimensional binary random variables, such that $0<E(X_{1})=\theta_{x}<\infty,\
0<Var(X_{1})=\theta_{x}(1-\theta_{x})<\infty$ and
$0<E(Y_{1})=\theta_{y}<\infty,\ 0<Var(Y_{1})=\theta_{y}(1-\theta_{y})<\infty$.
Further, let $S_{n}^{x}=\sum_{i=1}^{n}X_{i},\ S_{n}^{y}=\sum_{i=1}^{n}Y_{i}$
and let $M_{\delta}\equiv
M^{\delta}_{x}\equiv\inf\\{n>k_{x}^{\delta}:S_{n}^{x}>k_{x}^{\delta}\\}$
(which is denoted as $\tau(t)$ in [5]) be the corresponding stopping time.

######  Lemma 1.

Let $\eta^{2}:=Var_{\underset{{}^{\sim}}{\theta}}(\theta_{y}X-\theta_{x}Y)$,
since the joint distribution of $(X,Y)$ can be represented as in (3) of
Section $2$, we have by (4) that,

$\eta^{2}=\theta_{y}^{2}Var_{\theta_{x}}\left(X\right)+\theta_{x}^{2}Var_{\theta_{y}}\left(Y\right)-2\theta_{x}\theta_{y}Cov_{\underset{{}^{\sim}}{\theta}}\left(X,Y\right)=\theta_{x}\theta_{y}\left(\theta_{x}+\theta_{y}-2p_{11}\right)>0.$

Hence, with $M_{\delta}\equiv M^{\delta}_{x}$, and $\delta\to 0$,

$\frac{S^{y}_{M_{\delta}}}{k_{x}^{\delta}+1}\xrightarrow{P}\frac{\theta_{y}}{\theta_{x}},$

and the asymptotic distribution of $(S_{M_{\delta}}^{y},M_{\delta})^{\prime}$
is the normal ${\mathcal{N}}_{2}(\bm{\mu_{2}},\bm{V_{2}})$, where

$\bm{\mu_{2}}=\begin{pmatrix}[1.5]\frac{\theta_{y}}{\theta_{x}}\left(k^{\delta}_{x}+1\right)\\\
\frac{1}{\theta_{x}}\left(k^{\delta}_{x}+1\right)\end{pmatrix},\quad\bm{V_{2}}=\begin{pmatrix}[1.5]\frac{\theta_{y}\left(\theta_{x}+\theta_{y}-2p_{11}\right)\left(k^{\delta}_{x}+1\right)}{\theta_{x}^{2}}&\frac{\left(\theta_{y}-p_{11}\right)\left(k^{\delta}_{x}+1\right)}{\theta_{x}^{2}}\\\
\frac{\left(\theta_{y}-p_{11}\right)\left(k^{\delta}_{x}+1\right)}{\theta_{x}^{2}}&\frac{\left(1-\theta_{x}\right)\left(k^{\delta}_{x}+1\right)}{\theta_{x}^{2}}\end{pmatrix}.$

Similarly, when $M_{\delta}\equiv
M^{\delta}_{y}\equiv\inf\\{n>k_{y}^{\delta}:S_{n}^{y}>k_{y}^{\delta}\\}$, we
obtain as $\delta\to 0$,

$\frac{S^{x}_{M_{\delta}}}{k_{y}^{\delta}+1}\xrightarrow{P}\frac{\theta_{x}}{\theta_{y}},$

and $(S_{M_{\delta}}^{x},M_{\delta})^{\prime}$ is asymptotically normal
${\mathcal{N}}_{2}(\bm{\mu_{3}},\bm{V_{3}})$, where

$\bm{\mu_{3}}=\begin{pmatrix}[1.5]\frac{\theta_{x}}{\theta_{y}}\left(k^{\delta}_{y}+1\right)\\\
\frac{1}{\theta_{y}}\left(k^{\delta}_{y}+1\right)\end{pmatrix},\quad\bm{V_{3}}=\begin{pmatrix}[1.5]\frac{\theta_{x}\left(\theta_{x}+\theta_{y}-2p_{11}\right)\left(k^{\delta}_{y}+1\right)}{\theta_{y}^{2}}&\frac{\left(\theta_{x}-p_{11}\right)\left(k^{\delta}_{y}+1\right)}{\theta_{y}^{2}}\\\
\frac{\left(\theta_{x}-p_{11}\right)\left(k^{\delta}_{y}+1\right)}{\theta_{y}^{2}}&\frac{\left(1-\theta_{y}\right)\left(k^{\delta}_{y}+1\right)}{\theta_{y}^{2}}\end{pmatrix}.$

Now, by utilizing Lemma $1$, for small $\delta$ (as $\delta\to 0$), we obtain
that $\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(M_{\delta}=m)$, the
probability of the stopping time at $m$, for
$m=\overset{{}_{\sim}}{k}_{\delta}+1,\dots,N^{*}_{\delta}$, can be
approximated from the asymptotic bivariate normal distributions above, as

$\displaystyle\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{\delta}=m\right)\equiv\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\min\\{M_{x}^{\delta},M_{y}^{\delta}\\}=m\right)$
$\displaystyle=$
$\displaystyle\int_{-\infty}^{k^{\delta}_{y}+0.5}\int_{m-0.5}^{m+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{2}},\bm{V_{2}}\right)dwdu+\int_{-\infty}^{k^{\delta}_{x}+0.5}\int_{m-0.5}^{m+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{3}},\bm{V_{3}}\right)dwdu.$
(23)

In Appendix A.2 we provide all details leading to the approximation of
$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(M_{\delta}=m)$ in (5) above.
Immediately, by (5), we can approximate the values of
$ASN^{*}(\underset{{}^{\sim}}{\theta})$ in (11) and
$E_{\underset{{}^{\sim}}{\theta}}({M^{*}}^{2})$ in (16) as $\delta\to 0$. In a
similar manner, the power function
$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$ in (22) can
also be approximated, as $\delta\to 0$, as

$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})=\int_{-\infty}^{k^{\delta}_{y}+0.5}\int_{-\infty}^{N^{*}_{\delta}+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{2}},\bm{V_{2}}\right)dwdu+\int_{-\infty}^{k^{\delta}_{x}+0.5}\int_{-\infty}^{N^{*}_{\delta}+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{3}},\bm{V_{3}}\right)dwdu.$

######  Remark 1.

Note that as $\delta\to 0$, for our proposed test, we have

$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{x}^{\delta}<M_{y}^{\delta}\right)=\int_{-\infty}^{k^{\delta}_{y}+0.5}\int_{-\infty}^{N^{*}_{\delta}+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{2}},\bm{V_{2}}\right)dwdu,$

and

$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{y}^{\delta}<M_{x}^{\delta}\right)=\int_{-\infty}^{k^{\delta}_{x}+0.5}\int_{-\infty}^{N^{*}_{\delta}+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{3}},\bm{V_{3}}\right)dwdu.$

Accordingly, by (21) and (5) we can simplify the calculation of the expected
value of our post-test (post-detection) estimators,
$E_{\underset{{}^{\sim}}{\theta}}(\hat{\theta}_{x})$. Since we have, as
$\delta\to 0$,

$E_{\underset{{}^{\sim}}{\theta}}\left(\frac{S_{N^{*}_{\delta}}^{x}}{N^{*}_{\delta}}\mathbbm{1}\left[M_{\delta}>N^{*}_{\delta}\right]\right)=\int_{-\infty}^{k^{\delta}_{x}+0.5}\int_{-\infty}^{k^{\delta}_{y}+0.5}\frac{S_{N^{*}_{\delta}}^{x}}{N^{*}_{\delta}}\phi_{2}\left(u,w\mid\bm{\mu_{1}},\bm{V_{1}}\right)dwdu,$
(24)

and

$E_{\underset{{}^{\sim}}{\theta}}\left(\frac{S_{M^{*}_{\delta}}^{x}}{M^{*}_{\delta}}\mathbbm{1}\left[M_{\delta}\leq
N^{*}_{\delta}\right]\right)=\int_{-\infty}^{k^{\delta}_{y}+0.5}\int_{-\infty}^{N^{*}_{\delta}+0.5}\frac{k^{\delta}_{x}+1}{M_{\delta}}\phi_{2}\left(u,w\mid\bm{\mu_{2}},\bm{V_{2}}\right)dwdu\\\
+\int_{-\infty}^{k^{\delta}_{x}+0.5}\int_{-\infty}^{N^{*}_{\delta}+0.5}\frac{S_{M_{\delta}}^{x}}{M_{\delta}}\phi_{2}\left(u,w\mid\bm{\mu_{3}},\bm{V_{3}}\right)dwdu,$
(25)

we have a simplified expression of
$E_{\underset{{}^{\sim}}{\theta}}(\hat{\theta}_{x})$ by combining (24) and
(25) as $\delta\to 0$. And also for
$E_{\underset{{}^{\sim}}{\theta}}(\hat{\theta}_{y})$ as $\delta\to 0$,

$E_{\underset{{}^{\sim}}{\theta}}\left(\frac{S_{N^{*}_{\delta}}^{y}}{N^{*}_{\delta}}\mathbbm{1}\left[M_{\delta}>N^{*}_{\delta}\right]\right)=\int_{-\infty}^{k^{\delta}_{x}+0.5}\int_{-\infty}^{k^{\delta}_{y}+0.5}\frac{S_{N^{*}_{\delta}}^{y}}{N^{*}_{\delta}}\phi_{2}\left(u,w\mid\bm{\mu_{1}},\bm{V_{1}}\right)dwdu,$
(26)

and

$E_{\underset{{}^{\sim}}{\theta}}\left(\frac{S_{M^{*}_{\delta}}^{y}}{M^{*}_{\delta}}\mathbbm{1}\left[M_{\delta}\leq
N^{*}_{\delta}\right]\right)=\int_{-\infty}^{k^{\delta}_{x}+0.5}\int_{-\infty}^{N^{*}_{\delta}+0.5}\frac{k^{\delta}_{y}+1}{M_{\delta}}\phi_{2}\left(u,w\mid\bm{\mu_{2}},\bm{V_{2}}\right)dwdu\\\
+\int_{-\infty}^{k^{\delta}_{y}+0.5}\int_{-\infty}^{N^{*}_{\delta}+0.5}\frac{S_{M_{\delta}}^{y}}{M_{\delta}}\phi_{2}\left(u,w\mid\bm{\mu_{3}},\bm{V_{3}}\right)dwdu,$
(27)

we have the easier expression of
$E_{\underset{{}^{\sim}}{\theta}}(\hat{\theta}_{y})$ by combining (26) and
(27) as $\delta\to 0$.

Reflective of the effects of the ’stopping time’ $M^{*}_{\delta}$ on the
parameter estimates $\hat{\theta}_{x}$ (and $\hat{\theta}_{y}$), we do not
expect these estimators to generally be unbiased. However, to study the
question of the possible bias of the post-test (post-detection) estimator
(e.g. $\hat{\theta}_{x}$), we may consider the relative absolute bias between
the estimate and the actual value should be

$\frac{\lvert
E_{\underset{{}^{\sim}}{\theta}}\left(\hat{\theta}_{x}\right)-\theta_{x}\rvert}{\theta_{x}}\cdot
100\%,$

where $\lvert\cdot\rvert$ denotes the absolute value. Additionally, we denote
by, $\gamma={\theta_{x}}/{\theta_{y}}$, the relative risk of the two side
effects, $X$ and $Y$, and by $\gamma_{0}=\theta_{x}^{0}/\theta_{y}^{0}$, its
value under the null hypothesis.

In the following Table 4 and Figure 4, we present the maximal relative
absolute bias and its tendency between the estimated and the actual value,
which also illustrates the results of Theorem 2 below which indicates the
magnitude of the bias tends to $0$. In this illustration, we separate three
cases of the relative risk, $\gamma_{0}=0.5,1,2$. From Figure 4, we can
conclude that the relative risk does influence the tendency of the maximal
relative absolute bias. Since the computing time of the exact calculation
based on (18)-(20) substantially increases as the sample size increases, we
use the exact calculation only when $\delta>0.5$ and otherwise we use the
Monte Carlo approximation to evaluate (24) and (25). For instance, from Table
4, we can see that when $\delta=0.6$ and $\gamma_{0}=0.5$, the maximal
relative absolute bias is only $1.8528\%$, which is, in practical terms, very
small.

Table 4: The percent $(0-100)$ of the maximal relative absolute bias of $E_{\underset{{}^{\sim}}{\theta}}(\hat{\theta}_{x})$ from exact calculation and Monte Carlo approximation. We assume that (i) $\theta_{x}^{0}=0.05$, $\theta_{y}^{0}=0.1$; (ii) $\theta_{x}^{0}=\theta_{y}^{0}=0.05$; (iii) $\theta_{x}^{0}=0.1$, $\theta_{y}^{0}=0.05$. In all cases, $\alpha_{x}=\alpha_{y}=0.025$, $\beta_{x}=\beta_{y}=0.1$ and $\rho=0.1$. |  maximal relative absolute bias $(\%)$
---|---
scenario | $\delta=1$ | $\delta=0.8$ | $\delta=0.6$ | $\delta=0.4$ | $\delta=0.2$ | $\delta=0.1$
$\gamma_{0}=0.5$ | 3.8618 | 2.8326 | 1.8528 | 4.4194 | 2.5015 | 1.3984
$\gamma_{0}=1$ | 4.5492 | 3.2657 | 2.0921 | 4.4093 | 2.5868 | 1.4513
$\gamma_{0}=2$ | 4.1259 | 3.1066 | 2.0177 | 4.3121 | 2.4608 | 1.3853

(a) $\gamma_{0}=0.5$

(b) $\gamma_{0}=1$

(c) $\gamma_{0}=2$

Figure 4: The plots of the tendency of the maximal relative absolute bias
among various values of $\theta_{x}$.

######  Theorem 2.

Let $\hat{\theta}_{x}$ and $\hat{\theta}_{y}$ be the post-test (post-
detection) estimators (in (17)) of $\theta_{x}$ and $\theta_{y}$ based on the
$M^{*}_{\delta}=\min\left\\{N^{*}_{\delta},M_{\delta}\right\\}$ observations
obtained from the optimal $(\alpha,\beta)$ bivariate sequential test with
$N_{\delta}^{*}$, $k_{x}^{\delta}$ and $k_{y}^{\delta}$. Then as $\delta\to
0$, we have $\forall\ \underset{{}^{\sim}}{\theta}\in\Theta_{R}$,

$\hat{\theta}_{x}\xrightarrow{P}\theta_{x},\quad\text{and}\quad\hat{\theta}_{y}\xrightarrow{P}\theta_{y}.$

Proof We start with $\hat{\theta}_{x}$. For small $\delta$ (as $\delta\to 0$),
since we can express $\hat{\theta}_{x}$ by combining the following three terms
together, which is

$\hat{\theta}_{x}=\frac{S^{x}_{N^{*}_{\delta}}}{N^{*}_{\delta}}\mathbbm{1}\left[M_{\delta}>N^{*}_{\delta}\right]+\frac{k^{\delta}_{x}+1}{M_{\delta}}\mathbbm{1}\left[M_{\delta}\leq
N^{*}_{\delta}\text{ and
}M_{x}^{\delta}<M_{y}^{\delta}\right]+\frac{S^{x}_{M_{\delta}}}{M_{\delta}}\mathbbm{1}\left[M_{\delta}\leq
N^{*}_{\delta}\text{ and }M_{x}^{\delta}>M_{y}^{\delta}\right].$ (28)

The first term in (28) expresses the case that we would not reject the null
hypothesis, in which case $S^{x}_{N^{*}_{\delta}}$ is a binomial process with
$N^{*}_{\delta}$ and $\theta_{x}$. Hence, we have

$E_{\underset{{}^{\sim}}{\theta}}\left(\frac{S^{x}_{N^{*}_{\delta}}}{N^{*}_{\delta}}\right)=\theta_{x}\text{
and }\lim_{\delta\to
0}Var_{\underset{{}^{\sim}}{\theta}}\left(\frac{S^{x}_{N^{*}_{\delta}}}{N^{*}_{\delta}}\right)=\lim_{\delta\to
0}\frac{\theta_{x}\left(1-\theta_{x}\right)}{N^{*}_{\delta}}=0.$

which indicates that as $\delta\to 0$,
$\hat{\theta}_{x}\xrightarrow{P}\theta_{x}$ under the case that we would not
reject the null hypothesis.

The second term in (28) expresses the case that we would reject the null
hypothesis by stopping the process at the boundary $k_{x}^{\delta}+1$
corresponding to the side effect $X$ since $S_{M_{\delta}}^{x}>k^{\delta}_{x}$
and $S_{M_{\delta}}^{y}\leq k^{\delta}_{y}$. From the proof of Theorem $5$ in
[17], as $\delta\to 0$,

$\hat{\theta}_{x}\equiv\frac{k^{\delta}_{x}+1}{M_{\delta}}\xrightarrow{P}\theta_{x},$
(29)

which indicates that if we reject the null hypothesis upon stopping the
process at the boundary $k_{x}^{\delta}+1$,
$\hat{\theta}_{x}\xrightarrow{P}\theta_{x}$ as $\delta\to 0$.

The third term in (28) expresses the case that we reject the null hypothesis
by stopping the process at the boundary $k_{y}^{\delta}+1$ corresponding to
the side effect $Y$ since $S_{M_{\delta}}^{y}>k^{\delta}_{y}$ and
$S_{M_{\delta}}^{x}\leq k^{\delta}_{x}$. Since
$S_{M_{\delta}}^{y}>k^{\delta}_{y}$, by the result obtaining in Lemma 1, we
have as $\delta\to 0$,

$\frac{S_{M_{\delta}}^{x}}{k^{\delta}_{y}+1}\xrightarrow{P}\frac{\theta_{x}}{\theta_{y}}.$

Again, from the proof of Theorem $5$ in [17], we have as $\delta\to 0$,

$\frac{k^{\delta}_{y}+1}{M_{\delta}}\xrightarrow{P}\theta_{y}.$ (30)

Hence, we obtain

$\hat{\theta}_{x}\equiv\frac{S_{M_{\delta}}^{x}}{M_{\delta}}\xrightarrow{P}\theta_{x},\quad\text{as
}\delta\to 0.$

Accordingly, in all these three cases, we have as $\delta\to 0$,

$\hat{\theta}_{x}\xrightarrow{P}\theta_{x}.$

In a similar manner, we prove that as $\delta\to 0$, $\forall\
\underset{{}^{\sim}}{\theta}\in\Theta_{R}$,

$\hat{\theta}_{y}\xrightarrow{P}\theta_{y}.$

$\Box$

Once Theorem 2 is established, we now introduce one of our main results about
the joint asymptotic normality of the post-test (post-detection) estimator
$\underset{{}^{\sim}}{\hat{\theta}}\equiv(\hat{\theta}_{x},\hat{\theta}_{y})^{\prime}$.
This result enables us to construct an approximate $(1-\alpha)\cdot 100\%$
joint confidence interval of the real $\underset{{}^{\sim}}{\theta}$.

######  Theorem 3.

$\forall$ $\underset{{}^{\sim}}{\theta}\in\Theta_{R}$, as $\delta\rightarrow
0$, we have

$\bm{\mathcal{U}}_{\delta}:=\sqrt{M^{*}_{\delta}}\left(\underset{{}^{\sim}}{\hat{\theta}}-\underset{{}^{\sim}}{\theta}\right)\xrightarrow{D}\mathcal{N}_{2}\left(\bm{0},{\bm{\Sigma}}\right),$
(31)

where

${\bm{\Sigma}}=\begin{pmatrix}[1.5]\theta_{x}(1-\theta_{x})&p_{11}-\theta_{x}\theta_{y}\\\
p_{11}-\theta_{x}\theta_{y}&\theta_{y}(1-\theta_{y})\end{pmatrix}.$

Proof We have for any
$\underset{{}^{\sim}}{t}=\left(t_{1},t_{2}\right)^{\prime}\in\mathbb{R}^{2}$,

$\displaystyle\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\bm{\mathcal{U}}_{\delta}\leq\underset{{}^{\sim}}{t}\right)$
$\displaystyle=\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\bm{\mathcal{U}}_{\delta}\leq\underset{{}^{\sim}}{t}\mid
M_{\delta}>N^{*}_{\delta}\right)\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{\delta}>N^{*}_{\delta}\right)$
$\displaystyle+\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\bm{\mathcal{U}}_{\delta}\leq\underset{{}^{\sim}}{t}\mid
M_{\delta}\leq N^{*}_{\delta}\text{ and
}M_{x}^{\delta}<M_{y}^{\delta}\right)\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{\delta}\leq
N^{*}_{\delta}\text{ and }M_{x}^{\delta}<M_{y}^{\delta}\right)$
$\displaystyle+\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\bm{\mathcal{U}}_{\delta}\leq\underset{{}^{\sim}}{t}\mid
M_{\delta}\leq N^{*}_{\delta}\text{ and
}M_{x}^{\delta}>M_{y}^{\delta}\right)\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{\delta}\leq
N^{*}_{\delta}\text{ and }M_{x}^{\delta}>M_{y}^{\delta}\right)$
$\displaystyle:=B_{1}\cdot P_{1}+B_{2}\cdot P_{2}+B_{3}\cdot P_{3}.$

Hence,

$\lim_{\delta\rightarrow
0}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\bm{{\mathcal{U}}}_{\delta}\leq\underset{{}^{\sim}}{t}\right)=\lim_{\delta\rightarrow
0}B_{1}\cdot\lim_{\delta\rightarrow 0}P_{1}+\lim_{\delta\rightarrow
0}B_{2}\cdot\lim_{\delta\rightarrow 0}P_{2}+\lim_{\delta\rightarrow
0}B_{3}\cdot\lim_{\delta\rightarrow 0}P_{3}.$ (32)

Note that $P_{1}+P_{2}+P_{3}\to 1$, as $\delta\to 0$. The term $B_{1}$ in (32)
expresses the case that we would not reject the null hypothesis, which means
$M^{*}_{\delta}\equiv N^{*}_{\delta}$ observations. In which case, we have
$\underset{{}^{\sim}}{\hat{\theta}}=(\frac{S_{N^{*}_{\delta}}^{x}}{N^{*}_{\delta}},\frac{S_{N^{*}_{\delta}}^{y}}{N^{*}_{\delta}})^{\prime}$.
By the asymptotic distribution of
$(S_{N^{*}_{\delta}}^{x},S_{N^{*}_{\delta}}^{y})^{\prime}$ in (21), we
immediately obtain that

$\bm{\mathcal{U}}_{\delta}\equiv\sqrt{N^{*}_{\delta}}\left(\underset{{}^{\sim}}{\hat{\theta}}-\underset{{}^{\sim}}{\theta}\right)\xrightarrow{D}\mathcal{N}_{2}\left(\bm{0},{\bm{\Sigma}}\right),\quad\text{as}\quad\delta\to
0.$

Considering now the second term, $B_{2}$, in (32), we notice that we would
stop early at the boundary $S_{M_{\delta}}^{x}=k^{\delta}_{x}+1$ corresponding
to the side effect $X$, in which case, the stopping time $M_{\delta}\equiv
M_{x}^{\delta}$.

By utilizing the asymptotic bivariate normality of
$(S_{M_{\delta}}^{y},M_{\delta})^{\prime}$ in Lemma 1, we may derive the
distribution of

$\underset{{}^{\sim}}{\hat{\theta}}\equiv
g_{1}\left(S_{M_{\delta}}^{y},M_{\delta}\right):=\begin{pmatrix}[1.5]\frac{k^{\delta}_{x}+1}{M_{\delta}}\\\
\frac{S_{M_{\delta}}^{y}}{M_{\delta}}\end{pmatrix},$

by applying the standard delta method. Accordingly, we compute the
corresponding gradient vector of the function $g_{1}$ and evaluate it at
$\bm{\mu_{2}}$. Therefore, the asymptotic variance is

$\nabla g_{1}\bm{V_{2}}{\nabla
g_{1}}^{\prime}=\begin{pmatrix}[1.5]\frac{\theta_{x}^{2}(1-\theta_{x})}{k^{\delta}_{x}+1}&\frac{\theta_{x}\left(p_{11}-\theta_{x}\theta_{y}\right)}{k^{\delta}_{x}+1}\\\
\frac{\theta_{x}\left(p_{11}-\theta_{x}\theta_{y}\right)}{k^{\delta}_{x}+1}&\frac{\theta_{x}\theta_{y}(1-\theta_{y})}{k^{\delta}_{x}+1}\end{pmatrix}.$

Hence, we obtain that as $\delta\to 0$,

$\sqrt{\frac{k^{\delta}_{x}+1}{\theta_{x}}}\left(\underset{{}^{\sim}}{\hat{\theta}}-\underset{{}^{\sim}}{\theta}\right)\xrightarrow{D}\mathcal{N}_{2}\left(\bm{0},{\bm{\Sigma}}\right),\quad\text{
as }\quad{\delta}\rightarrow 0.$

By the result stated in (29), we conclude that in this case (when
$M_{\delta}\leq N^{*}_{\delta}$ and $M_{x}^{\delta}<M_{y}^{\delta}$),

$\bm{\mathcal{U}}_{\delta}\equiv\sqrt{M_{\delta}}\left(\underset{{}^{\sim}}{\hat{\theta}}-\underset{{}^{\sim}}{\theta}\right)\xrightarrow{D}\mathcal{N}_{2}\left(\bm{0},{\bm{\Sigma}}\right),\quad\text{
as }\quad{\delta}\rightarrow 0.$

To consider the third term $B_{3}$ in (32), we notice that we would stop early
at the boundary $S_{M_{\delta}}^{y}=k^{\delta}_{y}+1$ corresponding to the
side effect $Y$, in which case, the stopping time $M_{\delta}\equiv
M_{y}^{\delta}$.

By utilizing the asymptotic bivariate normality of
$(S_{M_{\delta}}^{x},M_{\delta})^{\prime}$ in Lemma 1, we may derive the
distribution of

$\underset{{}^{\sim}}{\hat{\theta}}\equiv
g_{2}\left(S_{M_{\delta}}^{x},M_{\delta}\right):=\begin{pmatrix}[1.5]\frac{S_{M_{\delta}}^{x}}{M_{\delta}}\\\
\frac{k^{\delta}_{y}+1}{M_{\delta}}\end{pmatrix},$

by applying the standard delta method. Accordingly, we compute the
corresponding gradient vector of the function $g_{2}$ and evaluate it at
$\bm{\mu_{3}}$. Therefore, the asymptotic variance is,

$\nabla g_{2}\bm{V_{3}}{\nabla
g_{2}}^{\prime}=\begin{pmatrix}[1.5]\frac{\theta_{x}\theta_{y}(1-\theta_{x})}{k^{\delta}_{y}+1}&\frac{\theta_{y}\left(p_{11}-\theta_{x}\theta_{y}\right)}{k^{\delta}_{y}+1}\\\
\frac{\theta_{y}\left(p_{11}-\theta_{x}\theta_{y}\right)}{k^{\delta}_{y}+1}&\frac{\theta_{y}^{2}(1-\theta_{y})}{k^{\delta}_{y}+1}\end{pmatrix}.$

Hence, we obtain that as $\delta\to 0$,

$\sqrt{\frac{k^{\delta}_{y}+1}{\theta_{y}}}\left(\underset{{}^{\sim}}{\hat{\theta}}-\underset{{}^{\sim}}{\theta}\right)\xrightarrow{D}\mathcal{N}_{2}\left(\bm{0},{\bm{\Sigma}}\right),\quad\text{
as }\quad{\delta}\rightarrow 0.$

By the result stated in (30), we conclude that in this case (when
$M_{\delta}\leq N^{*}_{\delta}$ and $M_{x}^{\delta}>M_{y}^{\delta}$),

$\bm{\mathcal{U}}_{\delta}\equiv\sqrt{M_{\delta}}\left(\underset{{}^{\sim}}{\hat{\theta}}-\underset{{}^{\sim}}{\theta}\right)\xrightarrow{D}\mathcal{N}_{2}\left(\bm{0},{\bm{\Sigma}}\right),\quad\text{
as }\quad{\delta}\rightarrow 0.$

Hence, upon combining the above together (and accounting of
$P_{1}+P_{2}+P_{3}\to 1$ as $\delta\to 0$), we obtain

$\displaystyle\lim_{\delta\rightarrow
0}\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left({\cal
U}_{\delta}\leq\underset{{}^{\sim}}{t}\right)$
$\displaystyle=\Phi_{2}\left(\underset{{}^{\sim}}{t}\mid\bm{0},{\bm{\Sigma}}\right)\cdot\lim_{\delta\rightarrow
0}P_{1}+\Phi_{2}\left(\underset{{}^{\sim}}{t}\mid\bm{0},{\bm{\Sigma}}\right)\cdot\lim_{\delta\rightarrow
0}P_{2}+\Phi_{2}\left(\underset{{}^{\sim}}{t}\mid\bm{0},{\bm{\Sigma}}\right)\cdot\lim_{\delta\rightarrow
0}P_{3}$
$\displaystyle=\Phi_{2}\left(\underset{{}^{\sim}}{t}\mid\bm{0},{\bm{\Sigma}}\right),$

where $\Phi_{2}(\underset{{}^{\sim}}{t}\mid\bm{0},{\bm{\Sigma}})$ denotes the
joint $cdf$ of the bivariate normal distribution
$\mathcal{N}_{2}(\bm{0},{\bm{\Sigma}})$. Accordingly, we have that

$\sqrt{M^{*}_{\delta}}\left(\underset{{}^{\sim}}{\hat{\theta}}-\underset{{}^{\sim}}{\theta}\right)\xrightarrow{D}\mathcal{N}_{2}\left(\bm{0},{\bm{\Sigma}}\right)\quad\text{
as }\quad{\delta}\rightarrow 0,$

as stated in (31). $\Box$

We close this section with a brief discussion of the (asymptotic) properties
of the estimated relative risk between the two side effects, $X$, and $Y$,
namely of $\gamma:=\theta_{x}/\theta_{y}$. Indeed, the following result is an
immediate consequence to the results stated in Theorem 3 and it deals with the
asymptotic distribution of
$\hat{\gamma}:={\hat{\theta}_{x}}/{\hat{\theta}_{y}}$.

######  Theorem 4.

Let $\hat{\gamma}$ be as above. Then $\forall\
\underset{{}^{\sim}}{\theta}\in\Theta_{R}$, as $\delta\rightarrow 0$, we have

$\sqrt{M^{*}_{\delta}}\left(\hat{\gamma}-\gamma\right)\xrightarrow{D}\mathcal{N}\left(0,\gamma\left(\frac{\gamma+1}{\theta_{y}}-\frac{2p_{11}}{\theta_{y}^{2}}\right)\right).$

Outline of the Proof Since $\hat{\gamma}\equiv
h(\hat{\theta}_{x},\hat{\theta}_{y})$, by applying the standard delta method
to the result in Theorem 3, we compute the corresponding gradient vector of
the function $h$ and evaluate it at $\underset{{}^{\sim}}{\theta}$. Therefore,
the asymptotic variance is

$\nabla h{\bm{\Sigma}}{\nabla
h}^{\prime}=\frac{\theta_{x}(\theta_{x}+\theta_{y}-2p_{11})}{\theta_{y}^{3}}\equiv\gamma\left(\frac{\gamma+1}{\theta_{y}}-\frac{2p_{11}}{\theta_{y}^{2}}\right).$

Hence, the stated result follows.

Specifically, when we stop at the boundary at $k^{\delta}_{x}+1$ corresponding
to the side effect $X$, by (29), we have as $\delta\to 0$,

$\sqrt{k^{\delta}_{x}+1}\left(\hat{\gamma}-\gamma\right)\xrightarrow{D}\mathcal{N}\left(0,\gamma^{2}\left(\gamma+1-\frac{2p_{11}}{\theta_{y}}\right)\right);$

when we stop at the boundary at $k^{\delta}_{y}+1$ corresponding to the side
effect $Y$, by (30), we have as $\delta\to 0$,

$\sqrt{k^{\delta}_{y}+1}\left(\hat{\gamma}-\gamma\right)\xrightarrow{D}\mathcal{N}\left(0,\gamma\left(\gamma+1-\frac{2p_{11}}{\theta_{y}}\right)\right).$

$\Box$

######  Remark 2.

Similarly, denote $\nu:=1/\gamma$. We may obtain the asymptotic normality of
$\hat{\nu}:=1/\hat{\gamma}$. That is, $\forall\
\underset{{}^{\sim}}{\theta}\in\Theta_{R}$, as $\delta\to 0$,

$\sqrt{M^{*}_{\delta}}\left({\hat{\nu}}-{\nu}\right)\equiv\sqrt{M^{*}_{\delta}}\left(\frac{1}{\hat{\gamma}}-\frac{1}{\gamma}\right)\xrightarrow{D}\mathcal{N}\left(0,\nu\left(\frac{\nu+1}{\theta_{x}}-\frac{2p_{11}}{\theta_{x}^{2}}\right)\right).$

##  6 Analysis of Some COVID-19 Side Effects Data

[16] provided the data on the side effects to COVID-19 vaccine which were
recorded among some health care workers in Nigeria. Their study accounted for
$117$ participants who received the COVID-19 vaccine. These vaccinated
participants reported on several side effects, if any. In the following Table
5, we provide the counts of some of the reported side effects, as fever,
muscle pain, dizziness, headache, etc.

|  | fever |  |
---|---|---|---|---
|  | No | Yes |
muscle pain | No | 63 | 11 |
| Yes | 18 | 25 | 43
|  |  | 36 | 117
(a)

|  | dizziness |  |
---|---|---|---|---
|  | No | Yes |
headache | No | 78 | 5 |
| Yes | 26 | 8 | 34
|  |  | 13 | 117
(b)

Table 5: $2\times 2$ table of example $1$ and example $2$

(a) example 1

(b) example 2

Figure 6: The plots of the $95\%$ simultaneous confidence intervals of example
$1$ and example $2$.

Example 1:

For this example, we focus attention on the side effect muscle pain, $X$, and
the side effect fever, $Y$. The results of these classifications are provided
as a $2\times 2$ table in Table 5 (a) above. For the post-test joint inference
of $(\theta_{x},\theta_{y})$, we have to estimate $\theta_{x}$, $\theta_{y}$
and $p_{11}$ in order to estimate $\bm{\Sigma}$ in (31). According to (3) and
(4), we calculate $\hat{\theta}_{x}=0.3675$, $\hat{\theta}_{y}=0.3077$, and by
direct calculation $\hat{p}_{11}=\frac{n_{11}}{M^{*}}=\frac{25}{117}=0.2137$
and the corresponding $\hat{\rho}=0.4521$. Further, the $95\%$ confidence
region is the ellipse shown in Figure 6 (a), which is centered at
$(0.3675,0.3077)$ and the half-lengths of the major and minor axes are
$0.1288$ and $0.0790$. By Theorem 2 and Theorem 3, we obtain the $95\%$
simultaneous confidence intervals are $\theta_{x}\in[0.2584,0.4766]$,
$\theta_{y}\in[0.2032,0.4121]$ and the $95\%$ Bonferroni confidence intervals
are $\theta_{x}\in[0.3008,0.4343]$, $\theta_{y}\in[0.2438,0.3716]$.

In the following, we assume two possible illustrative scenarios which could
have yielded the outcomes reported in table 5 (a) .

$\bullet$ Scenario I (upon rejection of $H_{0}$) We reconstruct the optimal
$(\alpha=0.05,\beta=0.09)$ sequential test of (1). To that end, we assume here
that $\theta_{x}^{0}=\theta_{y}^{0}=0.1$, then the pair of the test hypotheses
are:

$H_{0}:\theta_{x}\leq 0.1\text{ and }\theta_{y}\leq 0.1\ \ \text{against}\ \
H_{1}:\theta_{x}>0.1\text{ or }\theta_{y}>0.1.$

According to (7) and (8), assuming the termination occurring at $X$ boundary,
$M_{x}\equiv 117<N^{*}$, $M^{*}=\min\\{M_{x},N^{*}\\}=117$ with $k^{*}_{x}=42$
and $S_{M_{x}}=43$. With $\alpha_{x}=0.025$ and the probability of Type II
error, $\beta_{x}=0.09$, in particular, we assume $\theta_{x}^{1}=0.16$, we
obtain the maximal sample size $N^{*}_{x}=324$. Then for the same setup of
side effect $Y$, we obtain $N_{y}^{*}=324$, $k^{*}_{y}=42$. With no doubt that
$N^{*}_{x}=N^{*}_{y}$, we have $N^{*}=\min\\{N^{*}_{x},N^{*}_{y}\\}=324$.
Since the estimated $\hat{\rho}=0.4521$, we calculate the probabilities of
Type I error and Type II error in practice as

$\Pr\left(\text{Type~{}I~{}error}\right)=0.0561\quad\text{and}\quad\Pr\left(\text{Type~{}II~{}error}\right)=0.0208.$

Therefore, for this optimal $(\alpha\approx 0.05,\beta\approx 0.09)$
sequential test, with total sample size $N^{*}\equiv 324$ and the critical
value $k^{*}_{x}=k^{*}_{y}=42$, we would stop the study and reject the null
hypothesis.

$\bullet$ Scenario II (a non-rejection of $H_{0}$) We reconstruct the optimal
$(\alpha=0.05,\beta=0.092)$ sequential test of (1). But now we assume
$\theta_{x}^{0}=\theta_{y}^{0}=0.4$, then the pair of the test hypotheses are:

$H_{0}:\theta_{x}\leq 0.4\text{ and }\theta_{y}\leq 0.4\ \ \text{against}\ \
H_{1}:\theta_{x}>0.4\text{ or }\theta_{y}>0.4.$

According to (7) and (8), assuming $M=\min\\{M_{x},M_{y}\\}>N^{*}\equiv 117$,
we have $M^{*}=\min\\{M,N^{*}\\}=117$ and $S_{N^{*}}^{x}=43$,
$S_{N^{*}}^{y}=36$. Assume $X$ and $Y$ have the same setup, we have
$N^{*}=N_{x}^{*}=N_{y}^{*}=117$. With $\alpha_{x}=0.025$, we calculate the
critical value $k^{*}_{x}=57$. To achieve the probability of Type II error,
$\beta_{x}=0.092$, in particular, we assume $\theta_{x}^{1}=0.55$. Then for
the same setup of side effect $Y$, (that is,
$\theta_{y}^{0}=0.4,\theta_{y}^{1}=0.55,\alpha_{y}=0.025,\beta_{y}=0.092$), we
obtain $N^{*}_{y}=117$, $k^{*}_{y}=57$. Since the estimated
$\hat{\rho}=0.4521$, we calculate the probabilities of Type I error and Type
II error in practice as

$\Pr\left(\text{Type~{}I~{}error}\right)=0.0402\quad\text{and}\quad\Pr\left(\text{Type~{}II~{}error}\right)=0.0302.$

Therefore, for this optimal $(\alpha\approx 0.05,\beta\approx 0.092)$
sequential test, with total sample size $N^{*}\equiv 117$ and the critical
value $k^{*}_{x}=k^{*}_{y}=57$, we would not reject the null hypothesis that
the probabilities of participants exhibit muscle pain and exhibit fever both
are less than or equal to $0.4$.

Moreover, in this case, we may estimate the relative risk between the two side
effects by $\hat{\gamma}=\hat{\theta}_{x}/\hat{\theta}_{y}=1.1944$. Based on
Theorem 2 and Theorem 4, utilizing the above estimators, we obtain the $95\%$
confidence interval is $[0.8740,1.5149]$. Since this confidence interval
includes the value of $\gamma=1$, it indicates that, for $\alpha=0.05$, we
would not reject the null hypothesis, $H_{0}:\gamma=1$, that the relative risk
between muscle pain and fever is $1$ (i.e. that $\theta_{x}=\theta_{y}$).

Example 2:

For this example, we focus attention on the side effect headache, $X$, and the
side effect dizziness, $Y$. The results of these classifications are
summarized as a $2\times 2$ table in Table 5 (b). For the post-test joint
inference of $(\theta_{x},\theta_{y})$, we have to estimate $\theta_{x}$,
$\theta_{y}$ and $p_{11}$ in order to estimate $\bm{\Sigma}$ in (31).
According to (3) and (4), we calculate $\hat{\theta}_{x}=0.2906$,
$\hat{\theta}_{y}=0.1111$, and by direct calculation
$\hat{p}_{11}=\frac{n_{11}}{M^{*}}=\frac{8}{117}=0.0684$ and the corresponding
$\hat{\rho}=0.2529$. Further, the $95\%$ confidence region is the ellipse
shown in Figure 6 (b), which is centered at $(0.2906,0.1111)$ and the half-
lengths of the major and minor axes are $0.1055$ and $0.0670$. By Theorem $2$
and Theorem $3$, we obtain the $95\%$ simultaneous confidence intervals are
$\theta_{x}\in[0.1879,0.3933]$, $\theta_{y}\in[0.0400,0.1822]$ and the $95\%$
Bonferroni confidence intervals are $\theta_{x}\in[0.2278,0.3534]$,
$\theta_{y}\in[0.0676,0.1546]$.

Similarly, in the following, we assume two possible illustrative scenarios
which could have yielded the outcomes reported in Table 5 (b).

$\bullet$ Scenario I (upon rejection of $H_{0}$) We reconstruct the optimal
$(\alpha=0.05,\beta=0.09)$ sequential test of (1). To that end, we assume here
$\theta_{x}^{0}=\theta_{y}^{0}=0.1$, then the pair of the test hypotheses are:

$H_{0}:\theta_{x}\leq 0.1\text{ and }\theta_{y}\leq 0.1\ \ \text{against}\ \
H_{1}:\theta_{x}>0.1\text{ or }\theta_{y}>0.1.$

According to (7) and (8), assuming the termination occurring at $X$ boundary,
$M_{x}\equiv 117<N^{*}$, $M^{*}=\min\\{M_{x},N^{*}\\}=117$ with $k^{*}_{x}=33$
and $S_{M_{x}}=34$. With $\alpha_{x}=0.025$ and the probability of Type II
error, $\beta_{x}=0.09$, in particular, we assume $\theta_{x}^{1}=0.17$, we
obtain the maximal sample size $N^{*}_{x}=243$. Then for the same setup of
side effect $Y$, we obtain $N_{y}^{*}=243$, $k^{*}_{y}=33$. With no doubt that
$N^{*}_{x}=N^{*}_{y}$, we have $N^{*}=\min\\{N^{*}_{x},N^{*}_{y}\\}=243$.
Since the estimated $\hat{\rho}=0.2529$, we calculate the probabilities of
Type I error and Type II error in practice as

$\Pr\left(\text{Type~{}I~{}error}\right)=0.0472\quad\text{and}\quad\Pr\left(\text{Type~{}II~{}error}\right)=0.0167.$

Therefore, for this optimal $(\alpha\approx 0.05,\beta\approx 0.09)$
sequential test, with total sample size $N^{*}\equiv 243$ and the critical
value $k^{*}_{x}=k^{*}_{y}=33$, we would stop the study and reject the null
hypothesis.

$\bullet$ Scenario II (a non-rejection of $H_{0}$) We reconstruct the optimal
$(\alpha=0.05,\beta=0.11)$ sequential test of (1). But now we assume
$\theta_{x}^{0}=\theta_{y}^{0}=0.31$, then the pair of the test hypotheses
are:

$H_{0}:\theta_{x}\leq 0.31\text{ and }\theta_{y}\leq 0.31\ \ \text{against}\ \
H_{1}:\theta_{x}>0.31\text{ or }\theta_{y}>0.31.$

According to (7) and (8), assuming $M=\min\\{M_{x},M_{y}\\}>N^{*}\equiv 117$,
we have $M^{*}=\min\\{M,N^{*}\\}=117$ and $S_{N^{*}}^{x}=34$,
$S_{N^{*}}^{y}=13$. Assume $X$ and $Y$ have the same setup, we have
$N^{*}=N_{x}^{*}=N_{y}^{*}=117$. With $\alpha_{x}=0.025$, we calculate the
critical value $k^{*}_{x}=46$. To achieve the probability of Type II error,
$\beta_{x}=0.11$, in particular, we assume $\theta_{x}^{1}=0.45$. Then for the
same setup of side effect $Y$, (that is,
$\theta_{y}^{0}=0.31,\theta_{y}^{1}=0.45,\alpha_{y}=0.025,\beta_{y}=0.11$), we
obtain $N_{y}^{*}=117$, $k^{*}_{y}=46$. Since the estimated
$\hat{\rho}=0.2529$, we calculate the probabilities of Type I error and Type
II error in practice as

$\Pr\left(\text{Type~{}I~{}error}\right)=0.0394\quad\text{and}\quad\Pr\left(\text{Type~{}II~{}error}\right)=0.0288.$

Therefore, for this optimal $(\alpha\approx 0.05,\beta\approx 0.11)$
sequential test, with total sample size $N^{*}\equiv 117$ and the critical
value $k^{*}_{x}=k^{*}_{y}=46$, we would not reject the null hypothesis that
the probabilities of participants exhibit muscle pain and exhibit fever both
are less than or equal to $0.31$.

Moreover, similarly, we may estimate the relative risk
$\hat{\gamma}=\hat{\theta}_{x}/\hat{\theta}_{y}=2.6154$. Based on Theorem 2
and Theorem 4, utilizing the above estimators, we obtain the $95\%$ confidence
interval is $[1.2578,3.9730]$. Since the confidence interval exceeds
$\gamma=1$, it indicates that, for $\alpha=0.05$, we would reject the null
hypothesis, $H_{0}:\gamma=1$, that the relative risk between headache and
dizziness is greater than $1$. In this case, we would conclude that people are
more likely to exhibit dizziness than headache after they received COVID-19
vaccine.

##  7 Summary and Discussion

In this paper, we develop an $(\alpha,\beta)$-optimal sequential testing
procedure for an early detection of two potential side effects of certain
treatment. This sequential testing procedure does not require the
specification of the correlation, $\rho$, (if any) between the two potential
side effects nor any assumptions concerning it. Our procedure assures that the
actual probabilities of Type I and Type II errors would not exceed some
desired levels of $(\alpha,\beta)$ for all the possible values of $\rho$.

Since there is no assumption on the value of the correlation, we utilize the
(’negative’ version of the) multinomial distribution, to derive the exact
expression of the ASN and the variance of the ’stopping time’ $M^{*}$.
However, some tight bounds on the ASN are shown to hold following some simpler
calculations, once some general information on $\rho$ is available. For
instance, if these two side effects are independent (so that $\rho=0$), we
have a simplified version of the ASN available. Following basic analysis of
the properties of the stopping time, we focus on the post-detection estimators
of the model’s parameters $\theta_{x}$ and $\theta_{y}$. We derive the exact
formulas for calculating the expectation of the post-test (post-detection)
estimators and similarly outline the derivation needed for calculating the
corresponding variance.

To offset the computing time needed for the exact calculations, especially for
values of $(\theta_{x},\theta_{y})$ in close neighborhood of
$(\theta_{x}^{0},\theta_{y}^{0})$, the asymptotic properties of the final
sample size (i.e. the stopping time) are important to analyze. We derive the
joint (bivariate) asymptotic normality of
$(S_{M_{\delta}}^{z},M_{\delta})^{\prime}$, for $z=x,y$, which is the crucial
result for our subsequent analyses. Based on this asymptotic distribution, we
approximate the probability distribution of the stopping time at each possible
value in its support; a distribution that we then utilize to calculate the
ASN, and the expectation and the variance of the post-test (post-detection)
estimators, etc.

Moreover, the large sample consistency and the joint asymptotic normality of
the post-detection estimators, enable us to also construct the asymptotic
normality of the estimated relative risk, $\gamma$, of the two side effects.
In Section $6$, we presented two examples (based on real-life data) involving
’non-detection’ and ’detection’ situations of the side effects. In both
examples, we apply our sequential testing procedure and calculate the post-
test estimators, the corresponding joint confidence intervals, and also the
estimated relative risk and its confidence interval. These examples clearly
illustrate that our procedure performs well in the two different scenarios we
assumed (such assumptions can be defined by specialists). We note that the
nominal probabilities of Type I and Type II errors in these examples are less
than the desired $(\alpha,\beta)$, since the critical values used do not
utilize the correlation between these two side effects. The nominal values of
these error probabilities are conservative in any applied situation. We
conclude these two examples with the construction of a significance test of
hypothesis concerning the relative risk, $\gamma$, between the two side
effects we are interested in.

We point out that in some situations (e.g. $\theta_{x}\ll\theta_{y}$), the
probability of $\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(M_{x}>M_{y})$
mentioned in Remark $1$ is close to $1$. It indicates that once we stop our
observation process and reject the null hypothesis, we are likely to stop at
the boundary $k_{y}^{*}+1$ corresponding to the side effect $Y$. In such a
case, the one-dimensional optimal sequential test proposed in [17] is
sufficient to detect the potentially significant side effect (specifically,
side effect $Y$). Similar conclusion can be obtained for the case when
$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(M_{x}<M_{y})\to 1$.

Also note that if there are more than two potential side effects to account
for, one can separate these side effects into multiple pairs (decision should
involve input from the specialists), then applying our methods into each pair
to have further analysis results. To match the large sample size requirement,
our proposed test method can be sustainably applied to post-marketing
surveillance data.

In summary, we have demonstrated that our proposed sequential testing
procedure is particularly useful for an early detection of multiple side
effects especially in emergency situations as during the rapid deployment of
the COVID-19 vaccination campaign. Furthermore, the properties (asymptotic and
’finite-sample’) of the the post-test estimates of the unknown prevalence of
the potential side effects are very useful for any subsequent analysis.

##  8 Appendix

A.1 For side effect $X$, to construct the fixed-sample UMP test of

$H_{0}:\theta_{x}\leq\theta_{x}^{0}\ \ \text{against}\ \
H_{1}:\theta_{x}>\theta_{x}^{0},$

since the indicator for side effect $X$ is Bernoulli random variable, we have
the following properties:

$\tilde{\alpha}:=\operatorname{P}_{\theta_{x}^{0}}(S_{N_{x}}>k_{x})\leq\alpha,$
(33)

and

$\tilde{\beta}(\theta_{x}^{1})=\operatorname{P}_{\theta_{x}^{1}}(S_{N_{x}}\leq
k_{x})\leq\beta.$ (34)

Given corresponding $(\alpha,\beta,\theta_{x}^{0},\theta_{x}^{1})$,
$(\theta_{x}^{1}>\theta_{x}^{0})$, we may simultaneously ’solve’ equations
(33) and (34) for $N_{x}$ and $k_{x}$ to obtain the optimal ’sample size’,
$N^{*}_{x}\equiv N_{x}(\alpha,\beta,\theta_{x}^{0},\theta_{x}^{1})$ and a
corresponding ’critical test value’, $k^{*}_{x}$, by either an iterative
procedure utilizing (33) and (34) and the Binomial $pmf$ or by the standard
Normal approximations to the Binomial probabilities333See conditions in
$(1)-(2)$ of [8]. are given by,

$N^{*}_{x}=\left[\left(\frac{z_{\alpha}\sqrt{\theta_{x}^{0}(1-\theta_{x}^{0})}+z_{\beta}\sqrt{\theta_{x}^{1}(1-\theta_{x}^{1})}}{\theta_{x}^{1}-\theta_{x}^{0}}\right)^{2}\right],$

and,

$k^{*}_{x}=\left[N^{*}_{x}(z_{\alpha}\sqrt{\frac{\theta_{x}^{0}(1-\theta_{x}^{0})}{N^{*}_{x}}}+\theta_{x}^{0})-\frac{1}{2}\right],$

where $[x]$ is the nearest integer value to $x$ and $z_{p}:=-\Phi^{-1}(p)$,
$\forall\ p\in(0,1)$ where $\Phi$ denotes the standard Normal $cdf$.

Similarly, we can obtain the corresponding $N^{*}_{y}$ and $k^{*}_{y}$ for
side effect $Y$ by constructing the marginal fixed-sample UMP test for given
the values of $(\alpha,\beta,\theta_{y}^{0},\theta_{y}^{1})$ as well.

######  Lemma 2.

$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$, the power
function of optimal bivariate sequential test, it follows immediately from (9)
that, $\forall\ \underset{{}^{\sim}}{\theta}\in\Theta_{R}$, is monotonically
increasing with respect to $\theta_{x}$ and $\theta_{y}$.

Proof

$\displaystyle\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$
$\displaystyle=\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\operatorname{T^{*}_{seq}}\text{
reject
}H_{0}\right)=\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M\leq
N^{*}\right)=1-\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M>N^{*}\right)$
$\displaystyle=1-\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(S_{N^{*}}^{x}\leq
k^{*}_{x}\text{ and }S_{N^{*}}^{y}\leq k^{*}_{y}\right)$
$\displaystyle=1-\sum_{z=0}^{\overset{{}_{\sim}}{k}}\sum_{i=0}^{k^{*}_{x}-z}\sum_{j=0}^{k^{*}_{y}-z}\frac{N^{*}!}{z!i!j!(N^{*}-z-i-j)!}p_{11}^{z}p_{10}^{i}p_{01}^{j}p_{00}^{N^{*}-z-i-j}.$
(35)

To explore the tendency of the power function, we need to discuss separately
with respect to $\theta_{x}$ or $\theta_{y}$. Once we fixed
$\theta_{x}=\theta_{x}^{\prime}$, we may express the power function denoted as
$\Pi_{\operatorname{T^{*}_{seq}}}(\theta_{x}^{\prime},\theta_{y})\equiv\Pi^{\prime}(\theta_{y})$.
Note that in this case, the multinomial probabilities of
$\Pi^{\prime}(\theta_{y})$ in (8) would be reduced to a binomial
probabilities, a case which has been discussed in [17]. Accordingly, it has
been established there (see Theorem 1) that $\Pi^{\prime}(\theta_{y})$ is
monotonically increasing with respect to $\theta_{y}$. Similarly, when we fix
the value of $\theta_{y}$, we can establish the monotonicity of the power
function with respect to $\theta_{x}$ as well. Hence, we have the
monotonically increasing property of the power function
$\Pi_{\operatorname{T^{*}_{seq}}}(\underset{{}^{\sim}}{\theta})$ with respect
to its ordinates $\theta_{x}$ and $\theta_{y}$. $\Box$

######  Remark 3.

A key ingredient in the calculation of (11) was the probability mass function
of $M\equiv\min\\{M_{x},M_{y}\\}$. By utilizing the ’negative’ version of the
multinomial distribution in (2), with $p_{00}=1-\theta_{x}-\theta_{y}+p_{11}$,
$p_{10}=\theta_{x}-p_{11}$, $p_{01}=\theta_{y}-p_{11}$,
$p_{11}=Cov(X,Y)+\theta_{x}\theta_{y}$, we obtain by direct derivation that

$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M=m\right)\equiv\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\min\\{M_{x},M_{y}\\}=m\right)=A_{1}+A_{2}+A_{3}+A_{4}+A_{5},$
(36)

where

$A_{1}=\sum_{i=1}^{\min\\{k^{*}_{x}+1,k^{*}_{y}\\}}\sum_{j=0}^{\min\\{k^{*}_{y}-i,m-k^{*}_{x}-1\\}}A^{(1)}_{i,j},\quad\quad
A_{2}=\sum_{i=0}^{\overset{{}_{\sim}}{k}}\sum_{j=0}^{\min\\{k^{*}_{y}-i,m-k^{*}_{x}-1\\}}A^{(2)}_{i,j},$
$A_{3}=\sum_{i=1}^{\min\\{k^{*}_{x},k^{*}_{y}+1\\}}\sum_{j=0}^{\min\\{k^{*}_{x}-z,m-k^{*}_{y}-1\\}}A^{(3)}_{i,j},\quad\quad
A_{4}=\sum_{i=0}^{\overset{{}_{\sim}}{k}}\sum_{j=0}^{\min\\{k^{*}_{x}-i,m-k^{*}_{y}-1\\}}A^{(4)}_{i,j},$
$A_{5}=\sum_{i=1}^{\overset{{}_{\sim}}{k}+1}A^{(5)}_{i}\mathbbm{1}\left[m\geq
k^{*}_{x}+k^{*}_{y}+2-i\right],$

and

$A^{(1)}_{i,j}=\frac{(m-1)!}{(i-1)!(k^{*}_{x}+1-i)!j!(m-k^{*}_{x}-1-j)!}p_{11}^{i}p_{10}^{k^{*}_{x}+1-i}p_{01}^{j}p_{00}^{m-k^{*}_{x}-1-j},$
$A^{(2)}_{i,j}=\frac{(m-1)!}{i!(k^{*}_{x}-i)!j!(m-k^{*}_{x}-1-j)!}p_{11}^{i}p_{10}^{k^{*}_{x}+1-i}p_{01}^{j}p_{00}^{m-k^{*}_{x}-1-j},$
$A^{(3)}_{i,j}=\frac{(m-1)!}{(i-1)!j!(k^{*}_{y}+1-i)!(m-k^{*}_{y}-1-j)!}p_{11}^{i}p_{10}^{j}p_{01}^{k^{*}_{y}+1-i}p_{00}^{m-k^{*}_{y}-1-j},$
$A^{(4)}_{i,j}=\frac{(m-1)!}{i!j!(k^{*}_{y}-i)!(m-k^{*}_{y}-1-j)!}p_{11}^{i}p_{10}^{j}p_{01}^{k^{*}_{y}+1-i}p_{00}^{m-k^{*}_{y}-1-j},$
$A^{(5)}_{i}=\frac{(m-1)!}{(i-1)!(k^{*}_{x}+1-i)!(k^{*}_{y}+1-i)!(m-k^{*}_{y}-k^{*}_{x}-2+i)!}p_{11}^{i}p_{10}^{k^{*}_{x}+1-i}p_{01}^{k^{*}_{y}+1-i}p_{00}^{m-k^{*}_{y}-k^{*}_{x}-2+i}.$

Here, $A_{1}$ and $A_{2}$ indicate that the process terminates because of
exhibiting too many cases of side effect $X$, where $A_{1}$ represents the
last observation exhibited side effect $X$ only and $A_{2}$ represents the
last observation exhibited both side effects $X$ and $Y$; $A_{3}$ and $A_{4}$
indicate that the process terminates because of exhibiting too many cases of
side effect $Y$, where $A_{3}$ represents the last observation exhibited side
effect $Y$ only and $A_{4}$ represents the last observation exhibited both
side effects $X$ and $Y$; $A_{5}$ is the special case when the random walk on
the lattice stopped at the corner $(k^{*}_{x}+1,k^{*}_{y}+1)$.

######  Lemma 3.

Let $X$ be a random variable with nature numbers, then if $X$ is curtailed by
fixed size $n$,

$E\left(X\mathbbm{1}[X\leq n]\right)=\sum_{i=1}^{n}\operatorname{P}(X\geq
i)-n\operatorname{P}(X\geq n+1).$ (37)

Proof If $X$ is curtailed by fixed size $n$:

$\displaystyle E\left(X\mathbbm{1}[X\leq n]\right)=$
$\displaystyle\sum_{i=0}^{n}i\operatorname{P}(X=i)$ $\displaystyle=$
$\displaystyle\operatorname{P}(X=1)+2\cdot\operatorname{P}(X=2)+3\cdot\operatorname{P}(X=3)+\dots+n\cdot\operatorname{P}(X=n)$
$\displaystyle=$ $\displaystyle\operatorname{P}(X=1)$ $\displaystyle+$
$\displaystyle\operatorname{P}(X=2)+\operatorname{P}(X=2)$ $\displaystyle+$
$\displaystyle\operatorname{P}(X=3)+\operatorname{P}(X=3)+\operatorname{P}(X=3)$
$\displaystyle\vdots$ $\displaystyle+$
$\displaystyle\operatorname{P}(X=n)+\operatorname{P}(X=n)+\operatorname{P}(X=n)+\dots+\operatorname{P}(X=n)$
$\displaystyle=$ $\displaystyle\sum_{i=1}^{n}\operatorname{P}(X\geq
i)-n\operatorname{P}(X\geq n+1).$

$\Box$

######  Lemma 4.

$\forall\ \underset{{}^{\sim}}{\theta}\in\Theta_{R}$, if the correlation
between side effect $X$ and side effect $Y$ is positive ($\rho>0$), the
stopping time $M_{x}$ and the stopping time $M_{y}$ are positively associated;
if $\rho<0$, $M_{x}$ and $M_{y}$ are negatively associated; if $\rho=0$,
$M_{x}$ and $M_{y}$ are independent.

Proof Notice that
$Cov(X,Y)=\rho\sqrt{\theta_{x}(1-\theta_{x})\theta_{y}(1-\theta_{y})}$, if
$\rho>0$, $Cov(X,Y)>0$; if $\rho<0$, $Cov(X,Y)<0$. Since
${(X_{i},Y_{i})}_{i=1}^{n}$ are i.i.d. two-dimensional random variables, we
have $\forall\ n_{1},\ n_{2}\in\mathbbm{N^{+}}$,

$Cov\left(S_{n_{1}}^{x},S_{n_{2}}^{y}\right)=Cov\left(\sum_{i=1}^{n_{1}}X_{i},\sum_{j=1}^{n_{2}}Y_{j}\right)=\sum_{i=1}^{\min\\{n_{1},n_{2}\\}}Cov\left(X_{i},Y_{i}\right).$

Hence, if $\rho>0$, we have $Cov\left(S_{n_{1}}^{x},S_{n_{2}}^{y}\right)>0$;
if $\rho<0$, we have $Cov\left(S_{n_{1}}^{x},S_{n_{2}}^{y}\right)<0$. Since by
definitions of $M_{x}$ and $M_{y}$ in (6) can be represent as

$M_{x}=\sum_{n=0}^{\infty}\mathbbm{1}\left[S_{n}^{x}\leq
k^{*}_{x}\right]\quad\text{and}\quad
M_{y}=\sum_{n=0}^{\infty}\mathbbm{1}\left[S_{n}^{y}\leq k^{*}_{y}\right],$

where $\mathbbm{1}[S_{n}^{x}\leq k^{*}_{x}]$ and $\mathbbm{1}[S_{n}^{y}\leq
k^{*}_{y}]$ are nonincreasing functions of $S_{n}^{x}$ and $S_{n}^{y}$.

Hence ’$-M_{x}$’ and ’$-M_{y}$’ are nondecreasing functions of $S_{n}^{x}$ and
$S_{n}^{y}$ and since $Cov(-M_{x},-M_{y})=Cov(M_{x},M_{y})$, we have if
$\rho>0$, $M_{x}$ and $M_{y}$ are positively associated; if $\rho<0$, $M_{x}$
and $M_{y}$ are negatively associated (for properties of associated random
variables see [3] and [6]). On the other hand, since when $\rho=0$,
$cov(X,Y)=0$, we have $p_{11}=\theta_{x}\theta_{y}$. Therefore, $X$ and $Y$
are independent. Then since $M_{x}$ and $M_{y}$ are Borel-measurable functions
of $X$ and $Y$, we have $M_{x}$ and $M_{y}$ are independent. $\Box$

######  Lemma 5.

Let $X$ be a random variable with nature numbers, then

$E\left(X^{2}\right)=2\sum_{i=1}^{\infty}\left(i-\frac{1}{2}\right)\operatorname{P}\left(X\geq
i\right).$

If $X$ is curtailed by fixed size $n$, then

$E\left(X^{2}\mathbbm{1}\left[X\leq
n\right]\right)=2\sum_{i=1}^{n}\left(i-\frac{1}{2}\right)\operatorname{P}\left(X\geq
i\right)-n^{2}\operatorname{P}\left(X\geq n+1\right).$

Proof For $X$ is a random variable with nature numbers, since

$\displaystyle\frac{E\left(X^{2}\right)+E\left(X\right)}{2}=$ $\displaystyle
E\left(\frac{X(X+1)}{2}\right)=\sum_{i=1}^{\infty}\frac{i(i+1)}{2}\operatorname{P}\left(X=i\right)$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{\infty}\left(\sum_{j=1}^{i}j\right)\operatorname{P}\left(X=i\right)=\sum_{j=1}^{\infty}j\sum_{i=j}^{\infty}\operatorname{P}\left(X=i\right)$
$\displaystyle=$ $\displaystyle\sum_{j=1}^{\infty}j\operatorname{P}\left(X\geq
j\right).$

If $X$ is curtailed by fixed size $n$, since

$\displaystyle E\left(\frac{X(X+1)}{2}\mathbbm{1}\left[X\leq n\right]\right)=$
$\displaystyle\sum_{j=1}^{n}j\operatorname{P}\left(X\geq j\right)$
$\displaystyle=$ $\displaystyle\operatorname{P}\left(X=1\right)$
$\displaystyle+$
$\displaystyle\operatorname{P}\left(X=2\right)+2\operatorname{P}\left(X=2\right)$
$\displaystyle+$
$\displaystyle\operatorname{P}\left(X=3\right)+2\operatorname{P}\left(X=3\right)+3\operatorname{P}\left(X=3\right)$
$\displaystyle\vdots$ $\displaystyle+$
$\displaystyle\operatorname{P}\left(X=n\right)+2\operatorname{P}\left(X=n\right)+3\operatorname{P}\left(X=n\right)+\dots+n\operatorname{P}\left(X=n\right)$
$\displaystyle=$ $\displaystyle\sum_{i=1}^{n}i\operatorname{P}\left(X\geq
i\right)-\sum_{i=1}^{n}i\operatorname{P}\left(X\geq n+1\right)$
$\displaystyle=$ $\displaystyle\sum_{i=1}^{n}i\operatorname{P}\left(X\geq
i\right)-\frac{n(n+1)}{2}\operatorname{P}\left(X\geq n+1\right),$

and by (37) of Lemma $3$, we have

$\displaystyle E\left(X^{2}\mathbbm{1}\left[X\leq n\right]\right)=$
$\displaystyle 2E\left(\frac{X(X+1)}{2}\mathbbm{1}\left[X\leq
n\right]\right)-E\left(X\mathbbm{1}\left[X\leq n\right]\right)$
$\displaystyle=$ $\displaystyle
2\left[\sum_{i=1}^{n}i\operatorname{P}\left(X\geq
i\right)-\frac{n(n+1)}{2}\operatorname{P}\left(X\geq
n+1\right)\right]-\left[\sum_{i=1}^{n}\operatorname{P}\left(X\geq
i\right)-n\cdot\operatorname{P}\left(X\geq n+1\right)\right]$ $\displaystyle=$
$\displaystyle
2\sum_{i=1}^{n}\left(i-\frac{1}{2}\right)\operatorname{P}\left(X\geq
i\right)-n^{2}\operatorname{P}\left(X\geq n+1\right).$

$\Box$

A.2 Note that in the evaluation of
$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{\delta}=m\right)$ in
(36) of Remark $3$ above, involve the evaluation of

$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{\delta}=m\right)=A_{1}+A_{2}+A_{3}+A_{4}+A_{5}$

can be written as, by Lemma $1$, as $\delta\to 0$,

$\displaystyle\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{\delta}=m\right)\equiv\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\min\\{M_{x}^{\delta},M_{y}^{\delta}\\}=m\right)$
$\displaystyle=$
$\displaystyle\int_{-\infty}^{k^{\delta}_{y}+0.5}\int_{m-0.5}^{m+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{2}},\bm{V_{2}}\right)dwdu+\int_{-\infty}^{k^{\delta}_{x}+0.5}\int_{m-0.5}^{m+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{3}},\bm{V_{3}}\right)dwdu+A_{5}.$
(38)

Note that the first term in (8) indicates the probability when the last
observation exhibiting side effect $X$ causes the sum of observations that
exhibited side effect $X$ to attain the critical value $k^{\delta}_{x}+1$ but
$S_{M_{\delta}}^{y}<k^{\delta}_{y}+1$ (corresponding to $A_{1}$ and $A_{2}$ in
(36)). Also note the second term in (8) indicates the probability when the
last observation exhibiting side effect $Y$ causes the sum of observations
that exhibited side effect $Y$ to attain the critical value $k^{\delta}_{y}+1$
but $S_{M_{\delta}}^{x}<k^{\delta}_{x}+1$ (corresponding to $A_{3}$ and
$A_{4}$ in (36)).

However, whenever $\delta$ is small (as $\delta\to 0$), the value of $A_{5}$
in (8) which presents the probability of the last terminal observation
exhibiting both side effect $X$ and side effect $Y$ causes both critical
values $k^{\delta}_{x}+1$ and $k^{\delta}_{y}+1$ to be attained, which can be
approximated as

$A_{5}=p_{11}\cdot\int_{k^{\delta}_{x}-0.5}^{k^{\delta}_{x}+0.5}\int_{k^{\delta}_{y}-0.5}^{k^{\delta}_{y}+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{4}},\bm{V_{4}}\right)dwdu,$

where

$\bm{\mu_{4}}=\begin{pmatrix}[1.5]\left(m-1\right)\theta_{x}\\\
\left(m-1\right)\theta_{y}\end{pmatrix},\quad\bm{V_{4}}=\begin{pmatrix}[1.5]\left(m-1\right)\theta_{x}(1-\theta_{x})&\left(m-1\right)\rho\sqrt{\theta_{x}(1-\theta_{x})\theta_{y}(1-\theta_{y})}\\\
\left(m-1\right)\rho\sqrt{\theta_{x}(1-\theta_{x})\theta_{y}(1-\theta_{y})}&\left(m-1\right)\theta_{y}(1-\theta_{y}).\end{pmatrix}$

where $\phi_{2}(u,w\mid\bm{\mu_{4}},\bm{V_{4}})$ denotes the bivariate normal
distribution of $(S_{m-1}^{x},S_{m-1}^{y})^{\prime}$.

When we transform $(S_{m-1}^{x},S_{m-1}^{y})^{\prime}$ as the bivariate normal
with mean $(0,0)^{\prime}$ and variance-covariance matrix as
$\begin{pmatrix}1&\rho\\\ \rho&1\end{pmatrix}$, we obtain the difference
values between upper bounds and lower bounds of the integration of
$S_{m-1}^{x}$ and $S_{m-1}^{y}$ are, as $\delta\to 0$,

$\frac{k_{x}^{\delta}+0.5-\left(m-1\right)\theta_{x}}{\sqrt{\left(m-1\right)\theta_{x}\left(1-\theta_{x}\right)}}-\frac{k_{x}^{\delta}-0.5-\left(m-1\right)\theta_{x}}{\sqrt{\left(m-1\right)\theta_{x}\left(1-\theta_{x}\right)}}=\frac{1}{\sqrt{\left(m-1\right)\theta_{x}\left(1-\theta_{x}\right)}}\to
0,$

since $m\geq\overset{{}_{\sim}}{k}_{\delta}+1\to\infty$ and similarly,

$\frac{k_{y}^{\delta}+0.5-\left(m-1\right)\theta_{y}}{\sqrt{\left(m-1\right)\theta_{y}\left(1-\theta_{y}\right)}}-\frac{k_{y}^{\delta}-0.5-\left(m-1\right)\theta_{y}}{\sqrt{\left(m-1\right)\theta_{y}\left(1-\theta_{y}\right)}}=\frac{1}{\sqrt{\left(m-1\right)\theta_{y}\left(1-\theta_{y}\right)}}\to
0.$

It indicates that the value of $A_{5}$ is negligible as $\delta\to 0$.
Therefore, as $\delta\to 0$, we may use the approximated expression of
$\operatorname{P}_{\underset{{}^{\sim}}{\theta}}(M_{\delta}=m)$, that is

$\displaystyle\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(M_{\delta}=m\right)\equiv\operatorname{P}_{\underset{{}^{\sim}}{\theta}}\left(\min\\{M_{x}^{\delta},M_{y}^{\delta}\\}=m\right)$
$\displaystyle=$
$\displaystyle\int_{-\infty}^{k^{\delta}_{y}+0.5}\int_{m-0.5}^{m+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{2}},\bm{V_{2}}\right)dwdu+\int_{-\infty}^{k^{\delta}_{x}+0.5}\int_{m-0.5}^{m+0.5}\phi_{2}\left(u,w\mid\bm{\mu_{3}},\bm{V_{3}}\right)dwdu,$

to simplify the calculations.

## References

  * [1] HO Hartley and ER Fitch  “A chart for the incomplete beta-function and the cumulative binomial distribution”  In _Biometrika_ 38.3/4 JSTOR, 1951, pp. 423–426
  * [2] Paul R Rider  “The negative binomial distribution and the incomplete beta function”  In _The American Mathematical Monthly_ 69.4 JSTOR, 1962, pp. 302–304
  * [3] James D Esary, Frank Proschan and David W Walkup  “Association of random variables, with applications”  In _The Annals of Mathematical Statistics_ 38.5 Institute of mathematical statistics, 1967, pp. 1466–1474
  * [4] Michael Woodroofe  “Nonlinear renewal theory in sequential analysis”  SIAM, 1982
  * [5] Allan Gut and Svante Janson  “The limiting behaviour of certain stopped sums and some applications”  In _Scandinavian Journal of Statistics_ JSTOR, 1983, pp. 281–292
  * [6] Kumar Joag-Dev and Frank Proschan  “Negative association of random variables with applications”  In _The Annals of Statistics_ JSTOR, 1983, pp. 286–295
  * [7] Albert W Marshall and Ingram Olkin  “A family of bivariate distributions generated by the bivariate Bernoulli distribution”  In _Journal of the American Statistical Association_ 80.390 Taylor & Francis, 1985, pp. 332–338
  * [8] Martin Schader and Friedrich Schmid  “Two rules of thumb for the approximation of the binomial distribution by the normal distribution”  In _The American Statistician_ 43.1 Taylor & Francis, 1989, pp. 23–24
  * [9] Christopher Jennison and Bruce W Turnbull  “Group sequential tests for bivariate response: Interim analyses of clinical trials with both efficacy and safety endpoints”  In _Biometrics_ JSTOR, 1993, pp. 741–752
  * [10] John Bryant and Roger Day  “Incorporating toxicity considerations into the design of two-stage phase II clinical trials”  In _Biometrics_ JSTOR, 1995, pp. 1372–1383
  * [11] Mark R Conaway and Gina R Petroni  “Bivariate sequential designs for phase II trials”  In _Biometrics_ JSTOR, 1995, pp. 656–664
  * [12] Mark R Conway  “Importance sampling to evaluate bivariate probabilities with an application to sequential Phase II trials”  In _Journal of Statistical Computation and Simulation_ 51.2-4 Taylor & Francis, 1995, pp. 249–262
  * [13] Hua Jin  “Alternative designs of phase II trials considering response and toxicity”  In _Contemporary clinical trials_ 28.4 Elsevier, 2007, pp. 525–531
  * [14] Chia-Min Chen and Yunchan Chi  “Curtailed two-stage designs with two dependent binary endpoints”  In _Pharmaceutical Statistics_ 11.1 Wiley Online Library, 2012, pp. 57–62
  * [15] Huan Yin, Weizhen Wang and Zhongzhan Zhang  “On construction of single-arm two-stage designs with consideration of both response and toxicity”  In _Biometrical Journal_ 61.6 Wiley Online Library, 2019, pp. 1462–1476
  * [16] Oluwatosin Ruth Ilori et al.  “The acceptability and side effects of COVID-19 vaccine among health care workers in Nigeria: a cross-sectional study”  In _F1000Research_ 10 F1000 Research Limited London, UK, 2022, pp. 873
  * [17] Jiayue Wang and Ben Boukai  “Early Detection of Treatment’s Side Effect: A Sequential Approach”  In _Available at SSRN 4704084_ , 2024