# A Seamless Phase I/II Platform Design with a Time-To-Event Efficacy Endpoint
for Potential COVID-19 Therapies

Thomas Jaki1,2 and Helen Barnett1 and Andrew Titman1 and Pavel Mozgunov1

1 Department of Mathematics and Statistics, Lancaster University, Lancaster,
UK.
2 MRC Biostatistics Unit, University of Cambridge, Cambridge, UK

###### Abstract

In the search for effective treatments for COVID-19, initial emphasis has been
on re-purposed treatments. To maximise the chances of finding successful
treatments, novel treatments that have been developed for this disease in
particular, are needed. In this manuscript we describe and evaluate the
statistical design of the AGILE platform, an adaptive randomized seamless
Phase I/II trial platform that seeks to quickly establish a safe range of
doses and investigates treatments for potential efficacy using a Bayesian
sequential trial design. Both single agent and combination treatments are
considered. We find that the design can identify potential treatments that are
safe and efficacious reliably with small to moderate sample sizes.

Keywords: Adaptive Platform Trial, Dose-Escalation, COVID-19, Randomized,
Seamless, Time-to-improvement.

##  1 Introduction

The emergence of COVID-19 and the ensuing pandemic has led to a widespread,
frantic, search for treatments. Despite large uncertainty about the underlying
pathogen and the natural history of disease, trials must start rapidly to
identify treatments to save lives, but also so that effective treatments can
be used in the response to the outbreak. A consequence of this is that trials
in COVID-19 in the first few weeks and months of the outbreak have focused on
re-purposed treatments (Cao et al., 2020; Wang et al., 2020; Horby et al.,
2020b).

While recently, some success with using re-purposed treatments has been
achieved (Beigel et al., 2020; Horby et al., 2020a), it is crucial that
development of treatments specifically developed for COVID-19 is also
undertaken in order to maximise the chances finding therapies to successfully
treat patients. The crucial difference of trials investigating novel therapies
(in contrast to re-purposed treatments) is that the range of safe and likely
effective doses is unknown. Therefore, an efficient dose-finding design
identifying safe and active doses to be studied in larger trials is essential.
While there exist a number of dose-finding design for early-phase dose-finding
trials evaluating toxicity and efficacy simultaneously (see e.g. Wages and
Tait, 2015; Mozgunov and Jaki, 2019b, and references therein), many of them
consider a binary efficacy endpoint (with few recent extension to continious
endpoints, see e.g. Hirakawa, 2012; Mozgunov and Jaki, 2019a). Time-to-event
endpoints with censoring at 28-days have, however, previously been used as a
clinically meaningful measure in a number of COVID-19 trials (Cao et al.,
2020; Wang et al., 2020; Beigel et al., 2020) and the argument has been made
that they should be considered in COVID-19 trials (Dodd et al., 2020).

While the majority of Phase I dose-finding trials are non-randomised, it is
agreed that in later phases, the gold standard for evaluating novel treatments
are well conducted randomized controlled clinical trials. At the same time, in
the light of the uncertainty about the adverse events caused by COVID-19, it
is essential to conduct a randomised dose-finding trial to ensure that the
risk of adverse events is correctly attributed to the drug under study rather
than to the disease itself. Moreover, it has been argued that adaptive designs
(e.g. Pallmann et al., 2018) are particularly suitable during a pandemic, also
in the light of the uncertainty about a novel disease (Stallard et al., 2020).
Therefore, a randomised adaptive dose-finding design evaluating both toxicity
and time-to-event efficacy would allow to answer the research questions of
interest in novel therapies for treating COVID-19.

It is also recognised that there is a number of novel therapies that have a
potential to be efficient in fighting the COVID-19. Therefore, it is crucial
to have a structure in place that would allow rapid enrolment of novel
therapies to ensure rapid decision-making, and, importantly, would allow for
effecient use of information between the studies, i.e. utilising the data from
the control treatment across different compounds. This can be achieved via a
platform trial (Meyer et al., 2020).

In this paper, we describe and evaluate the design developed and implemented
for the AGILE platform (Griffiths et al., 2020), an adaptive randomized
seamless Phase I/II dose-finding trial platform that seeks to quickly
establish a safe range of doses and investigates treatments for potential
efficacy using a Bayesian sequential trial design. The proposed design uses
statistical models to improve a decision making, and further efficiency is
gained by sharing control group patients between concurrent controls across
different candidates in the platform. We also extend the design for the trials
studying dual-agent combination of treatments.

The rest of manuscript proceeds as follows. Section 2 describes the platform
for single treatments while its performance is evaluated in simulations in
Section 3. The design for dual-agent combinations is proposed in Section 4 and
subsequently evaluated in Section 5. We conclude with a discussion (Section
6).

##  2 Single agent design

###  2.1 Setting

Consider a randomized controlled dose-escalation clinical trial in which $m$
increasing doses $d_{1}<d_{2}<\ldots<d_{m}$ of a single experimental treatment
are studied. Let ${d}_{0}=0$ be a dose of zero of the treatment, which is
subsequently referred to as the control arm (or, simply, control). The
inclusion of the control arm is motivated by the still emerging nature of the
symptoms associated with COVID-19 and the desire to avoid labelling potential
treatments as unsafe due to misclassifying non-treatment related symptoms. A
binary outcome of a random variable $Y$, $y=0$ is observed if no dose-limiting
event (DLE) is observed within $t_{safe}>0$ days after randomisation, and
$y=1$, otherwise. Let $p_{j}$ be the probability for a patient to experience a
DLE if given dose ${d}_{j}$. It is assumed that the risk of a DLE is a non-
decreasing function of dose, $p_{0}\leq p_{1}\leq\ldots\leq p_{m}$ and prior
information for the DLE probability of the control arm, $p_{0}$, is available.

As it is expected that the control arm is associated with a non-negative
(unknown) risk of DLE (or symptoms of the disease that can not be
distinguished from DLEs), the primary goal of the dose-escalation is
formulated in terms of the additional risk of a dose limiting event (ADLE)
defined in terms of the expected difference in DLE risk between the doses of
the agent and the control. Specifically, we therefore seek to identify the
dose that corresponds to an ADLE risk of $\gamma=0.20$ which equates to
finding the dose $d_{j^{\star}}$ such that

$j^{\star}=\arg\min_{j=0,\ldots,m}|\left(p_{j}-p_{0}\right)-\gamma|.$

###  2.2 Bayesian Dose-Escalation Model

The following randomized Bayesian dose-escalation design that builds on the
proposal by Mozgunov et al. (2019) is used. Assume that the DLE probability
has the functional form

$\psi(\tilde{d}_{j},\theta_{1},\theta_{2})=\frac{\exp(\theta_{1}+\theta_{2}\tilde{d}_{j})}{1+\exp(\theta_{1}+\theta_{2}\tilde{d}_{j})}$
(1)

where $\theta_{1}$ and $\theta_{2}$ are unknown parameters, and
$\tilde{d}_{j}$ is a standardized dose level (also referred to as a skeleton)
corresponding to dose $j$, which are constucted given the information about
the prior DLE toxicities (details are given below). This model choice was
found to result in good statistical properties in terms of the target dose
identification in a randomized dose-finding trial (Mozgunov et al., 2019). We
require that the standardized dose level corresponding to control is equal to
$\tilde{d}_{0}=0$. This will guarantee that a sequential update of the slope
parameter $\theta_{1}$ will not contribute to the DLE probability estimation
on the control arm yet all data are used for its estimation (Mozgunov et al.,
2019).

Denote the prior distribution of the vector $\theta=(\theta_{1},\theta_{2})$
by $f_{0}(.)$. To construct the standardized levels, $\tilde{d}_{j}$, we
represent them in terms of prior estimates of the DLE probabilities
$\hat{p}_{j}^{(0)}$ associated with doses ${d}_{j}$ $j=0,\ldots,m$

$\tilde{d}_{j}=\frac{{\rm
logit}(\hat{p}_{j}^{(0)})-\hat{\theta}_{1}^{(0)}}{\hat{\theta}_{2}^{(0)}}$ (2)

where $\hat{\theta}_{1}^{(0)}$, $\hat{\theta}_{2}^{(0)}$ are prior point
estimates of the model parameters, and ${\rm logit}(x)=\log\frac{x}{1-x}$ is
the logit transformation of $x$. To satisfy $\tilde{d}_{0}=0$, the prior needs
to be chosen such that ${\rm
logit}(\hat{p}_{0}^{(0)})=\hat{\theta}_{1}^{(0)}$.

Assume that $n$ patients have already been assigned to doses
$\tilde{d}(1),\ldots,\tilde{d}(n)$ and binary responses
$\mathbb{Y}_{n}=[y_{1},\ldots,y_{n}]^{\rm T}$ were observed, respectively. The
models updates the posterior distribution of $\theta$ using Bayes’ Theorem

$f_{n}(\theta)=\frac{f_{n-1}(\theta)\phi(\tilde{d}(n),y_{n},\theta)}{\int_{\mathbb{R}^{h}}f_{n-1}(u)\phi(\tilde{d}(n),y_{n},u){\rm
d}u}=\frac{f_{0}(\theta)\prod_{i=1}^{n}\phi(\tilde{d}(i),y_{i},\theta)}{\int_{\mathbb{R}^{h}}f_{0}(u)\prod_{i=1}^{n}\phi(\tilde{d}(i),y_{i},u){\rm
d}u}$ (3)

where
$\phi(d(n),y_{n},\theta)=\psi(d(n),\theta)^{y_{n}}(1-\psi(d(n),\theta))^{1-y_{n}}.$
This posterior distribution is then used to make the escalation/de-escalation
decision (Neuenschwander et al., 2008). Specifically, the first set of safe
doses is defined as the doses $d_{j}$ for which

$\mathbb{P}\left(p_{j}-p_{0}\geq\gamma+2\delta\right)<c_{\rm overdose}$ (4)

where $\gamma$ is the target ADLE risk, $\delta$ is the width of the interval
of DLE risk which we consider acceptable, $c_{\rm overdose}$ is the threshold
controlling overdosing, and the probability is found with respect to the
updated posterior distribution. Amongst the safe doses, the dose which
maximises

$\mathbb{P}\left(p_{j}-p_{0}\in[\gamma-\delta,\gamma+\delta]\right)$ (5)

is selected as the target dose.

###  2.3 Efficacy Design

####  2.3.1 Baeysian Efficacy Model

For the efficacy endpoint, a Cox proportional hazards model is assumed where
the hazard of recovery at time $t$ is given by $h(t;z)=h_{0}(t;z)\psi^{z}$
where $z=1$ corresponds to a treatment and $z=0$ to control. Initially the
cohort of patients who have graduated from the dose-escalation stage are
followed up for a total of $t_{eff}$ days. Based on their outcomes, a decision
is made to either stop for futility, stop for efficacy or recruit a further
cohort of patients. To improve power, controls recruited from other candidate
treatments or other doses of the same treatment are also used within the
evaluation, but this is restricted to using only the most recent $n_{c}$ such
controls to mitigate the risk of bias due to population drift.

A Bayesian criterion is adopted for the stopping rule at each stage. In line
with Bayesian thinking, we set the stopping rules to be the same for each
stage $k$. Specifically denoting all data up to stage $k$ on dose $j$ by
$\mathcal{D}^{(j)}_{k}$ and for a given desirable treatment effect,
$\psi^{*}$,

  * •

evaluation is stopped for efficacy if
$\pi^{(j)}_{E|k}\coloneqq\mathbb{P}(\psi=\psi^{*}\mid\mathcal{D}^{(j)}_{k})>u$,

  * •

evaluation is stopped for futility if $\pi^{(j)}_{E|k}<l$, or

  * •

an additional cohort of patients is recruited, otherwise.

At the final stage $k=K$, efficacy for dose $j$ is established if
$\pi^{(j)}_{E|K}>u$.

A point prior of the form
$\mathbb{P}(\psi=1)=1-\pi^{(j)}_{E},~{}\mathbb{P}(\psi=\psi^{*})=\pi^{(j)}_{E}$
is assumed for $\psi$. Here $\pi^{(j)}_{E}$ represents the degree of optimism
or scepticism towards the likely efficacy of dose $d_{j}$.

An advantage of the point prior is that obtaining the posterior probability
$\pi^{(j)}_{E|k}$ is computationally straightforward allowing comprehensive
evaluations of the design via simulations. The posterior probability under
this model is

$\pi^{(j)}_{E|k}=\frac{\pi^{(j)}_{E}L(\psi^{*}\mid\mathcal{D}^{(j)}_{k})}{\pi^{(j)}_{E}L(\psi^{*}\mid\mathcal{D}^{(j)}_{k})+(1-\pi^{(j)}_{E})L(1\mid\mathcal{D}^{(j)}_{k})},$
(6)

where $L(u\mid\mathcal{D}^{(j)}_{k})$ is the Cox partial likelihood with
respect to the data for dose $j$ up to the $k$th period evaluated at a hazard
ratio of $u$. While the Cox partial likelihood is often not considered
compatible with a Bayesian analysis since it does not use the full information
of the data, Bayesian justifications of its use are available (Kalbfleisch,
1978; Sinha et al., 2003).

####  2.3.2 Setting the boundaries

To set the boundaries, $(l,u)$, a large number of trajectories of
$L(1\mid\mathcal{D}^{(j)}_{k})$ and $L(\psi^{*}\mid\mathcal{D}^{(j)}_{k})$ can
be simulated under both the null and alternative hypothesis, where in all
cases the simulation continues until the maximum period $K$. An assumption
regarding the proportion of patients recovered by time $t_{eff}$ under the
null is needed, in addition to the hazard ratio between treatments. In the
absence of censoring due to drop-out, the results will be otherwise invariant
to the precise survival distribution assumed since the Cox partial likelihood
only uses the order of events.

The effect of varying $(l,u)$ can then be explored by converting the pairs of
likelihoods into a posterior probability and imposing the boundary stopping
rules. For any given set of boundaries, the type I error, power, expected
number of patients under the null, expected number of patients under the
alternative and probabilities of stopping for futility or efficacy at each
stage, can be approximated. The boundaries can then be set to optimize some
criterion, for instance sum of expected sample sizes under the null and
alternative, subject to some constraints, for instance controlling type I
error, keeping power above some level or limiting the chance of early stopping
under the alternative.

The inclusion of historic controls will increase both the power and type I
error of any procedure (e.g. Schmidli et al., 2014). To ensure type I error is
controlled, the boundaries are set assuming the maximum, $n_{c}$, previous
controls are available, with the consequence that the type I error will be
lower for the first treatment evaluated in the platform. This will also mean
the power will be lower for the first few evaluated treatments. However, given
the dose finding design, it is anticipated that the first evaluations will be
of less importance, as for safety reasons evaluations tend to start at sub-
optimal low doses.

###  2.4 Overall Design

Throughout the study, patients are allocated in cohorts of size
$c=c_{1}+c_{2}$ where $c_{1}$ is the number of patients in the cohort assigned
to an active dose and $c_{2}$ is the number of patients in cohort assigned to
the control arm, $\tilde{d}_{0}$. Below is an outline of the overall procedure
made up by both safety and efficacy evaluation:

Safety-evaluation

  1. 1.

The first cohort of $c_{1}+c_{2}$ patients is assigned to the first dose and
to the control arm, respectively.

  2. 2.

After $t_{safe}$ days, short-term DLE outcomes are collected and the posterior
distribution of the parameters is updated using Equation 3.

  3. 3.

The set of safe doses is found using Equation 4.

     * •

If no doses are safe, the trial is stopped for safety;

     * •

if only the current dose is safe, the next cohort of $c_{1}+c_{2}$ patients is
assigned to the current dose and to the control arm, respectively;

     * •

otherwise, the next cohort of $c_{1}+c_{2}$ patients is assigned to the
adjacent, safe dose level for which the probability (5) is maximized and
control arm respectively.

  4. 4.

Once efficacy information is available for two cohorts on a safe dose, that
dose is graduated to the efficacy evaluation.

Efficacy-evaluation

  1. 1.

If a dose $d_{j}$ is deemed safe, the efficacy outcome is observed up to day,
$t_{eff}$.

  2. 2.

The posterior probability, $\pi^{(j)}_{E|k}$, following equation 6 is then
computed where $k$ corresponds to the number of times this dose has been
evaluated for efficacy.

     * •

If $\pi^{(j)}_{E|k}<l$ evaluation of dose $d_{j}$ is stopped for futility;

     * •

if $\pi^{(j)}_{E|k}>u$ evaluation is stopped for efficacy and the
corresponding candidate and dose $d_{j}$ recommended for further testing;

     * •

otherwise if $k<K$, an additional cohort of $c_{1}+c_{2}$ patients is
recruited on the current dose and control arm, respectively.

The evaluation of a dose continues until the maximum number of patients $N$ on
a dose has been reached unless it is stopped for efficacy, futility or safety
before. Once all doses are stopped, the evaluation of this candidate stops.
The overall design of the study is depicted in Figure 1.

Figure 1: Illustration of the AGILE platform design.

##  3 Evaluations of Proposed Design

###  3.1 Setting

We will now evaluate safety and efficacy across the study together in a
simulation study and evaluate the impact of shared control data that are
gradually accumulated over the course of the trial.

We consider the setting where there are three active doses ($m=3$) of a single
agent and the control arm. As before, the DLE endpoint is binary and
corresponds to either experiencing or not experiencing a DLE by time
$t_{safe}=7$ and the efficacy endpoint is time-to-improvement defined as a 1
category improvement on the 10-point WHO scale (Marshall et al., 2020) over 28
days, $t_{eff}=28$. To generate the efficacy outcome a Weibull distribution is
used. For the control group the rate and shape parameters were 0.085 and
0.797, respectively, resulting in a median recovery time of 14 days and a
recovery rate of 70% within 28 days. The rate parameter for the candidate dose
have been adapted to match the scenarios described below. The binary safety
and time-to-event efficacy responses are assumed to be highly correlated and
are generated via a procedure described by Mozgunov et al. (2018) using a
correlation coefficient $\rho=0.8$.

The maximum total intake per dose level is 72 patients assigned to each dose
level and control, which equates to a maximum total sample size of 216. In
line with the real study, cohort sizes of $c_{1}=4$ (assigned to active doses)
and $c_{2}=2$ (assigned to the control) are used. Moreover, we also report
sensitivity analysis to assess the impact of various cohort size.

The objective of the trial is to find all safe efficacious doses to be
graduated into a larger Phase II or Phase III clinical trial. The target ADLE
risk is $\gamma=0.20$, and the dose is considered safe if the ADLE risk is
less than or equal to $0.30$ ($\delta=0.05$) and is efficacious if it
corresponds to hazard ratio of at least 1.75. The stopping boundaries for
efficacy have been found such that the type I error in each pairwise
comparison is controlled at 10% one-sided, while the power for each individual
dose is 80%. The resulting boundary values are $l=0.224$ and $u=0.839$ (see
Section 3.3.2 for more details).

###  3.2 Scenarios

As the trial aims to study novel compounds which have yet to be explored with
respect to their mechanism of action in COVID-19 patients, it is crucial that
the design has good operating characteristics under a variety of dose-DLE and
dose-efficacy scenarios. Therefore, we consider 5 dose-efficacy scenarios
ranging from no doses corresponding to a change in time-to-improvement within
28-days, to all doses resulting in a clinically significant reduction; and 5
dose-DLE scenarios ranging from all doses being safe to all doses being very
unsafe. We then consider all combinations of these scenarios, resulting in 25
scenarios explored in total. The five dose-DLE and dose-efficacy scenarios for
each for ($d_{0}$, $d_{1}$, $d_{2}$, $d_{3}$) are presented in Table 1.

Table 1: Safety and efficacy scenarios for ($d_{0}$, $d_{1}$, $d_{2}$,
$d_{3}$).

Safety (Probability of DLE, $p_{j}$) Efficacy (Hazard Ratio) Scenario 0 (0.10,
0.12, 0.13, 0.15) (1.00, 1.00, 1.00, 1.00) Scenario 1 (0.10, 0.12, 0.15, 0.30)
(1.00, 1.00, 1.75, 1.75) Scenario 2 (0.10, 0.15, 0.30, 0.45) (1.00, 1.50,
1.75, 1.75) Scenario 3 (0.10, 0.30, 0.45, 0.60) (1.00, 1.50, 1.75, 2.00)
Scenario 4 (0.10, 0.45, 0.60, 0.60) (1.00, 1.75, 2.00, 2.00)

We will refer to the scenario with dose-DLE relationship $x$ and dose-efficacy
relationship $y$ as “Scenario $x$–$y$”. Each dose under the combination of DLE
and efficacy scenario is classified as incorrect, undesirable, acceptable or
desirable. If a treatment is unsafe or has a hazard ratio of 1 then it is
classed as incorrect. If it is safe, then a hazard ratio of 1.25 is
undesirable, 1.5 is acceptable and at least 1.75 is desirable.

For all 25 scenarios, a sensitivity analysis is conducted on varying values of
$c_{1}$ and $c_{2}$ in order to assess the effect of both altering the
allocation ratio between control and active doses, and the total cohort size,
$c$. We also study the implications of not sharing controls across doses. A
total of six settings are considered, ($c_{1}=2$ & $c_{2}=1$ and $c_{1}=2$ &
$c_{2}=2$ for only the settings where controls are shared across doses, with
$c_{1}=4$ & $c_{2}=2$ and $c_{1}=3$ & $c_{2}=3$ for the settings both where
controls are shared across doses and where they are not). The maximum number
of cohorts per dose varies with cohort size in order to maintain the constant
maximum total intake per dose level of 72 across that dose level and control.

Software in the form of R code used to produce the presented results is
available on GitHub (https://github.com/dose-finding/covid19-agile).

###  3.3 Choice of design parameters

####  3.3.1 Safety Model

The proposed design requires the prior and design parameters for both safety
and efficacy parts to be pre-specified in advance of the conduct of the trial.
The procedure of how these parameters were chosen are given below.

The prior parameters for the safety model were obtained via a calibration
procedure (Lee and Cheung, 2009) over a number of safety scenarios (not taking
into account efficacy). We use safety scenarios 1–3 in Table 1 that correspond
to the target dose being $d_{3}$, $d_{2}$, and $d_{1}$, respectively, thus
covering various locations of the target dose on the dosage grid.

The following prior distribution for the vector of safety model parameters
$\theta$ was assumed:

$(\theta_{1},\log(\theta_{2}))\sim\mathcal{N}(\mu,\Sigma)$

where $\mu=\left(\mu_{1},\mu_{2}\right)^{\rm T}$ is the vector of means and

$\Sigma=\left[\begin{array}[]{ccc}\sigma_{1}&\sigma_{12}\\\
\sigma_{12}&\sigma_{2}\\\ \end{array}\right].$

Given the link between the prior toxicity on the control and the intercept
parameter $\theta_{1}$ as implied by Equation 2, $\mu_{1}={\rm
logit}(\hat{p}_{0}^{(0)})$ where the prior DLE probability at the control,
$\hat{p}_{0}^{(0)}$. Following discussions with the clinical team, the DLE
risk on control was set to $\hat{p}_{0}^{(0)}=0.10$. To reduce the
computational complexity of the calibration, the covariance between the model
parameters was assumed to be $\sigma_{12}=0.$ The rest of the parameters were
chosen by conducting simulations using various combinations of values of the
parameters on the grid,
$\mu_{2}=\\{-0.5,-0.25,0,0.25,0.5\\},\sigma_{1}=\\{1.2,1.3,1.4,1.5,1.6\\},\sigma_{2}=\\{0.15,0.25,0.35,0.45,0.55\\}.$

Furthermore, to define the standardized doses $\tilde{d}_{j}$ in Equation 2,
the prior toxicity probability at each dose should be assumed. As there is no
reliable information of the DLE rates in the COVID-19 population, the skeleton
was also calibrated. The grid of values for the prior toxicity risks are
chosen in terms of the difference in the probability of DLE between the
neighbouring doses. Specifically,
$\hat{p}_{j}^{(0)}=\hat{p}_{0}^{(0)}+\nu\times j,j=1,2,3$ where $\nu$ is the
difference in the toxicity probabilities between doses, which are then used to
find the skeleton using Equation 2. The grid of values of $\nu$ was included
as one of the parameters for the calibration
$\nu=\\{0.075,0.100,0.125,0.150\\}$. Below, we fix $c_{\rm overdose}=0.25$
that was previously found to result in good safety properties of the design
for the two-parameter logistic model (Bailey et al., 2009).

The calibration was performed as follows. For each combination of parameters
of $\nu,\mu_{2},\sigma_{1},\sigma_{2}$ on the specified of grid, 500
simulations were run under each of the three considered scenarios monitoring
the proportion of target dose selections. Then, the selected values of the
parameters are those that maximized the geometric mean (taken across scenarios
for the same combination of values of the parameters) of the proportion of the
target dose selection. This resulted in using
$\nu=0.125,\mu_{2}=-0.25,\sigma_{1}=1.40,\sigma_{2}=0.35$ for the further
design evaluation.

####  3.3.2 Efficacy Model

The efficacy stopping boundaries for a particular setting were taken as the
pair $(l,u)$ that maximizes

$\mbox{Criterion}=\mbox{Power}-\lambda\\{\mathbb{E}(N_{0})+\mathbb{E}(N_{1})\\},$
(7)

subject to the constraint that the type I error is less than or equal to 10%,
where $\mathbb{E}(N_{0})$ and $\mathbb{E}(N_{1})$ are the expected sample
sizes (across both active dose and control arms) under the null and
alternative hypothesis, respectively and $\lambda$ is a tuning parameter.

A value of $\lambda=1/320$ was chosen to allow a power of 80% to be achieved
when $\psi=1.75$ for the main settings where $n_{c}=30$ and $c_{1}:c_{2}=4:2$.
For comparability, the same $\lambda$ was used in the sensitivity analysis
where the size of cohorts and/or the use of past controls was varied. In these
other settings power is lower than 80%, ranging from 63% to 79%. The
boundaries for the scenarios are shown in Table 2.

Throughout, the point prior is taken to be $\pi^{(j)}_{E}=\frac{1}{2}$ for all
doses $j$. Note that if the boundaries are chosen based on (7), the operating
characteristics of the efficacy design are not affected by the choice of point
prior. Specifically, since (6) can be re-expressed in terms of the posterior
log-odds as

$\log\left(\frac{\pi^{(j)}_{E|k}}{1-\pi^{(j)}_{E|k}}\right)=\log\left(\frac{\pi^{(j)}_{E}}{1-\pi^{(j)}_{E}}\right)+\log\left(\frac{L(\psi^{*}|D^{(j)}_{k})}{L(1|D^{(j)}_{k})}\right),$

changing the value of $\pi^{(j)}_{E}$ merely has the effect of translating the
posterior log-odds by a constant. Hence the boundaries $(\tilde{l},\tilde{u})$
under the alternative prior will satisfy the relationship
$(\mbox{logit}(\tilde{l}),\mbox{logit}(\tilde{u}))=(\mbox{logit}(l)+\xi,\mbox{logit}(u)+\xi)$
where $\xi$ is the log-odds ratio between the new and old prior odds of
efficacy.

Table 2: Boundaries for settings in the sensitivity analysis of cohort sizes.
The settings have the same maximum sample sizes and common criterion to trade-
off power and average sample size.

Cohort $c_{1}$ $c_{2}$ $n_{c}$ $l$ $u$ Power Criterion 4 2 30 0.224 0.839
0.800 0.630 3 3 30 0.268 0.841 0.747 0.573 2 1 30 0.192 0.858 0.794 0.634 2 2
30 0.227 0.858 0.744 0.566 4 2 0 0.317 0.815 0.634 0.438 3 3 0 0.271 0.821
0.691 0.484

###  3.4 Results

Detailed results for the setting with cohort size of 4+2 and inclusion of up
to thirty control patients are presented in Table 3. A comparison of results
for the varying cohort sizes are illustrated in Figures 2 and 3. Across
settings, the type I error rate, that is the percentage of simulations in
Scenario 0–0 where any dose is recommended, ranges from 11% to 14%.

Table 3: Percentage of 10,000 simulations where each dose is recommended for
($d_{1}$, $d_{2}$, $d_{3}$) for $c_{1}=4$ & $c_{2}=2$, with controls shared
across doses. Desirable doses are highlighted in bold and Acceptable doses are
highlighted in italics. Note that these may sum to more than 100% for each
scenario as more than one dose can be recommended simultaneously.

Efficacy Scenario 0 1 2 3 4 Safety 0 (0.4, 2.3, 8.9) (0.4, 24.0, 67.2) (2.9,
24.2, 67.5) (2.9, 24.6, 80.6) (5.1, 29.7, 80.9) Scenario 1 (1.1, 5.9, 8.0)
(1.2, 47.4, 52.9) (8.6, 47.6, 53.6) (9.1, 47.3, 62.1) (14.0, 56.3, 61.6) 2
(5.5, 7.9, 2.4) (5.7, 48.4, 9.4) (31.7, 49.9, 9.8) (31.4, 49.3, 10.9) (45.0,
57.2, 10.9) 3 (8.4, 2.1, 0.1) (8.2, 7.3, 0.2) (40.4, 7.8, 0.2) (40.4, 7.9,
0.2) (55.0, 8.3, 0.2) 4 (5.0, 0.1, 0.0) (4.9, 0.4, 0.0) (17.6, 0.3, 0.0)
(17.6, 0.3, 0.0) (23.1, 0.4, 0.0)

Figure 2: Percentage of simulations that recommend all desirable doses (left)
and the percentage of simulations that recommend any desirable dose (right)
for different cohort size and composition and with and without sharing control
group data. Note that only 13 out of 25 efficacy/safety scenarios contain a
desirable dose. Figure 3: Average total sample size across simulations for all
scenarios.

Figure 2 shows the percentage of simulations where all desirable doses are
recommended (left) and where any desirable doses are recommended (right). For
the baseline setting of $c_{1}=4$ and $c_{2}=2$, the mean percentage of
simulations recommending all desirable doses is 28.2%, whilst 71.0% of
simulations recommended any desirable dose. As expected, safety/efficacy
scenarios where only one dose is desirable have higher percentages of
simulations recommending all desirable doses, whilst lower percentages of
simulations recommending any desirable dose. The highest power for
recommending any desirable dose is observed in the safety/efficacy scenarios
where all active doses are safe and efficacious, reaching 89% for Scenario
1–4.

The sensitivity analysis across the varying cohort settings shows there is a
very small difference in performance. The ordering of performance across
safety/efficacy scenarios is identical with only a small numerical difference.
However, we can see that not sharing controls across doses decreases the
performance. Power for recommending any desirable dose increases for
increasing cohort size $c$, on average increasing 6% for allocation 1:2 ($c=3$
and $c=6$) and 2% for allocation ratio 1:1 ($c=4$ and $c=6$). Again
unsurprisingly, a lower power is observed when controls are not shared across
doses, and within this a higher power when the allocation ratio between
control and active dose within cohort is equal. The highest power is achieved
for cohorts sizes $c_{1}=4$ and $c_{2}=2$ with controls being shared across
doses.

Figure 3 illustrates the average total sample size across the scenarios and
settings. On average about 65 patients are required in the setting used in the
trial with the total sample size exceeding 150 in only 1% of simulations
across scenarios. The scenarios with the smallest sample sizes are those where
all doses are unsafe and the trial is therefore stopped early for safety. In
such cases, it takes 30 patients across 6 weeks on average to reach the
conclusion of stopping early for safety. The scenarios with larger sample
sizes are those where all doses are safe and most are acceptable or only just
desirable (i.e. not the case where the hazard ratio is 2.00), as in these
cases more doses are taken to the efficacy part and more patients are required
to detect the smaller difference in hazard ratios.

It can be seen across settings that the larger the total cohort size, the
larger the total sample size. This also corresponds to the higher power
settings. When controls are shared across doses, altering the
control:treatment allocation from 1:1 to 1:2 decreases the total average
sample size, whereas when controls are not shared across doses this increases.
This also links to the relationship with the power in the corresponding
settings; a higher power is achieved for equal allocations when controls are
not shared. The sharing of controls means more patients can be allocated to
the novel treatment and hence more can be learnt about it. In all cases, the
average total sample size is below 90.

Table 3 gives more detail into which doses are recommended across simulations.
For example, in efficacy scenario 2 where the lowest dose is acceptable and
the higher two are desirable in terms of efficacy. In safety scenario 0 where
all doses are safe, the highest dose is chosen most often. In safety scenario
2 where only the highest dose is unsafe, the middle dose is chosen most often,
although less often than the highest desirable dose in scenario 0. In safety
scenario 4 where all doses are unsafe, the lowest does is chosen only 17.6% of
the time. It is clear that desirable doses are recommended most often, with
incorrect and undesirable doses rarely recommended. This gives insight that
the procedure is successful in identifying desirable doses of a single agent.

##  4 Extension to combination treatments

###  4.1 Setting

Consider now a randomized controlled dose-escalation dual-agent clinical trial
studying the combinations of $J$ doses $d_{1}<d_{2},\ldots,d_{J}$ of the first
compound (referred to as agent $A$) and of $L$ doses
$s_{1}<s_{2}<\ldots<s_{L}$ of the second compound (referred to as agent $B$).
As before, let $d_{0}=s_{0}$ be a zero dose of each compound, respectively,
correspondingly to the control treatment, and denote the combination of dose
$d_{j}$ of $A$ and dose $s_{l}$ of $B$ by $(d_{j},s_{l})$. Within each agent
(assuming the second agent is fixed), it is known prior to the trial that the
risk of a DLE monotonically increases with the dose. The objective of the
trial is then to study the safety of the combinations and to establish the
maximum tolerated combination, the combination corresponding to the target
ADLE over the control arm of 20%. Denoting the probability of DLE at
combination with doses $d_{j}$ and $s_{l}$ by $p(d_{j},s_{l})$, the
probability of each agent given individually by $p(d_{j})$ and $p(s_{l})$, and
the probability of DLE at the control by $p(d_{0},s_{0})=p_{0}$, the aim is to
find the combination $(d_{j^{\star}},s_{l^{\star}})$ minimising

$|\left(p(d_{j},s_{l})-p_{0}\right)-\gamma|.$

The fundamental difference to the single-agent setting introduced in Section
2.1 is that one cannot order all of the combinations of the compounds with
respect to the monotonically increasing risk of DLE despite the monotonicity
assumption being satisfied within each compound. For example, comparing
combination $(d_{1},s_{2})$ and $(d_{2},s_{1})$, the dose of one agent is
increased and another is decreased, and it is unknown prior to the trial which
of these effects prevails in the overall DLE risk associated with the
combinations. Consequently, the model-based design for monotherapies in
Section 2.2 cannot be used and an alternative dose-finding. Below, the
extension of the single-agent model-based design is suggested.

###  4.2 Dual-Agent Bayesian Dose-Escalation Model

For the considered randomized dual-agent combination setting, under the
assumption of independence of the compounds, the probability of a DLE
associated with combination $(d_{j},s_{l})$ can be written as

$p_{0}(d_{j},s_{l})=1-(1-p(d_{j}))(1-p(s_{l})).$ (8)

To allow for the interaction of the compounds in terms of the probability of a
DLE, we use the proposed model by Neuenschwander et al. (2015)

$odds(p(d_{j},s_{l}))=odds(p_{0}(d_{j},s_{l}))\times\exp\left(\eta\times\tilde{d}_{j}\times\tilde{s}_{l}\right),$

where $odds(p)=\frac{p}{1-p}$ is the odds transformation of the probability
$p$, and ${\eta}$ is the interaction coefficient, positive values of which
correspond to synergistic DLE risk, zero corresponds to additive effect
without interaction, and negative values correspond to antagonistic risk of
DLE, and $\tilde{d}_{j}$, $\tilde{s}_{l}$ are standardized dose levels
corresponding to the same DLE probability as the doses ${d}_{j},s_{l}$.

Note that $p(\cdot)$ in Equation (8) is the probability of a DLE associated
with one compound given as monotherapy as in Section 2.2. Therefore, we adopt
the 2-parameter logistic model given in Equation (1) for each agent
separately. Specifically, let

$p(d_{j})=\psi(\tilde{d}_{j},\theta_{1},\theta_{21})$ (9)

and

$p(s_{l})=\psi(\tilde{s}_{l},\theta_{1},\theta_{22})$ (10)

where
$\theta=(\theta_{1},\log\left(\theta_{21}\right),\log\left(\theta_{22}\right),\eta)$
are the unknown parameters with a Normal prior distribution,
$\theta\sim\mathcal{N}(\mu,\Sigma)$ where
$\mu=\left(\mu_{1},\mu_{21},\mu_{22},\mu_{\eta}\right)^{\rm T}$ is the vector
of means and

$\Sigma=\left[\begin{array}[]{cccc}\sigma_{1}&\sigma_{1,21}&\sigma_{1,22}&0\\\
\sigma_{1,21}&\sigma_{21}&0&0\\\ \sigma_{1,22}&0&\sigma_{22}&0\\\
0&0&0&\sigma_{\eta}\\\ \end{array}\right]$

As before, we require the standardized dose level corresponding to the control
treatment to be equal to $\tilde{d}_{0}=\tilde{s}_{0}=0$, so that the
intercept parameter of the 2-parameter model (8) relates to the probability of
DLE on the control only. Therefore, both single-agent models employ the same
intercept parameter $\theta_{1}$ as it corresponds to the probability of a DLE
at the same control treatment. Consequently, for small to moderate values of
probability of DLE, the intercept parameter in the single-agent model,
$\theta_{1}$ approximately equals the logit inverse-logit transformation of
the half of the probability of DLE on the control treatment subject ${\rm
logit}\left(\frac{p(d_{0},s_{0})}{2}\right)={\theta_{1}}$. This is used to
construct the standardized dose levels $\tilde{d}_{j}$, $\tilde{s}_{l}$ using
the prior means of the parameters $\theta$ and the prior probabilities of a
DLE at each combination similarly to the construction in Equation (2).

Parameters of the vector $\theta$ are the unknown quantities that define the
combination-DLE relationship. As in the single-agent design in Section 2.2,
the posterior distribution of these are sequentially updated using the data
collected during the trial using Bayes’ theorem. Specifically, denote the
joint prior distributions of vector $\theta$ by $f_{0}(.)$. Assume that $n$
patients have received the combinations $(d1),s(1))\ldots,(d(n),s(n))$ and
binary responses $\mathbb{Y}_{n}=[y_{1},\ldots,y_{n}]^{\rm T}$ were observed,
respectively. The models update the posterior distribution of $\theta$ as

$f_{n}(\theta)=\frac{f_{n-1}(\theta)\phi((d(n),s(n)),y_{n},\theta)}{\int_{\mathbb{R}^{4}}f_{n-1}(u)\phi((d(n),s(n),y_{n},u){\rm
d}u}=\frac{f_{0}(\theta)\prod_{i=1}^{n}\phi((d(i),s(i),y_{i},\theta)}{\int_{\mathbb{R}^{4}}f_{0}(u)\prod_{i=1}^{n}\phi((d(i),s(i),y_{i},u){\rm
d}u}$ (11)

where

$\phi((d(n),s(n),y_{n},\theta)=p(d(n),s(n),\theta)^{y_{n}}(1-p(d(n),s(n),\theta))^{1-y_{n}}.$

This posterior distribution is then used to make the escalation and de-
escalation decision during the trials as proposed below.

###  4.3 Dual-Agent Dose-Escalation Design

The above combination-DLE model is then used in the design in Section 2.4 in
place of the single-agent model. As the efficacy part of the dose-escalation
design proposed for monotherapies considered each dose individually, the
efficacy part of the combination study remains the same. Once the combination
of the compounds is established to be safe, it is graduated into the efficacy
part following the single-agent proposal and the same decision rules for
dropping for futility and safety.

##  5 Evaluation of combination treatment design

###  5.1 Scenarios

In order to evaluate the dual-agent design, we conduct a simulation study
comprising of scenarios with two dose levels of agent $A$ ($d_{1}$ & $d_{2}$)
and three dose levels of agent $B$ ($s_{1}$, $s_{2}$ & $s_{3}$). We consider
four dose-DLE scenarios ranging from a situation where all combinations are
safe to a case where all are unsafe. Four dose-efficacy scenarios ranging from
no efficacious combination to a steep monotonic within agent relationship are
considered yielding 16 safety-efficacy scenarios. The dual agent scenarios are
presented in Table 4 with definitions of incorrect, undesirable, acceptable
and desirable dose combinations remaining as they were previously defined for
single agent doses in Section 3.2. Here we fix cohort sizes to the previous
baseline setting of $c_{1}=4$ (assigned to active dose combinations) and
$c_{2}=2$ (assigned to the control) and allow controls to be shared across
dose combinations. The maximum number of patients per dose combination is 72
as before.

Table 4: Safety and efficacy scenarios for dual-agent combinations $A$ and
$B$. It is assumed the control arm remains with probability of DLE 0.10 and
hazard ratio 1.00.

Safety Efficacy (Probability of DLE, $p_{j}$) (Hazard Ratio) $d_{1}$ $d_{2}$
$d_{1}$ $d_{2}$ Scenario 0 $s_{1}$ 0.10 0.12 1.00 1.00 $s_{2}$ 0.13 0.15 1.00
1.00 $s_{3}$ 0.15 0.18 1.00 1.00 Scenario 1 $s_{1}$ 0.10 0.12 1.00 1.25
$s_{2}$ 0.25 0.30 1.25 1.50 $s_{3}$ 0.50 0.55 1.50 1.75 Scenario 2 $s_{1}$
0.15 0.30 1.00 1.50 $s_{2}$ 0.25 0.35 1.25 1.75 $s_{3}$ 0.30 0.45 1.50 2.00
Scenario 3 $s_{1}$ 0.40 0.45 1.00 1.50 $s_{2}$ 0.45 0.50 1.50 1.75 $s_{3}$
0.50 0.55 1.75 1.75

###  5.2 Parameter Calibration

To define the parameters of the combination model, a calibration procedure
similar to the procedure described in Section 3.3.1 was applied. Safety
scenarios 0–3 in Table 4 that correspond to different steepness of the
combination-toxicity relationship and different locations of the target
combination were used. We then choose the hyperparameters for the prior
distribution of the parameters of the model,
$\theta=(\theta_{1},\log\left(\theta_{21}\right),\log\left(\theta_{22}\right),\eta)$,
given in Equation (9) and Equation (10).

Given the link between the prior toxicity on the control arm, and the
intercept parameter $\theta_{1}$ that is common for both single-agent
parameter models, we set $\mu_{1}={\rm
logit}\left(\hat{p}_{0}^{(0)})/2\right)$ where $\hat{p}_{0}^{(0)}=0.10$ is the
DLE risk on the control, as before. To reduce the computational complexity,
the covariance between all model parameters was assumed to be
$\sigma_{1,21}=\sigma_{1,22}=0$, and the variance of the slope parameters in
each single-agent model is the same $\sigma_{21}=\sigma_{22}.$ Finally, we set
the mean of the distribution of the interaction parameter $\mu_{\eta}=0$ to
reflect that either synergistic and antogonistic effects are possible. The
rest of the parameters were chosen by simulations using combinations of values
on the grid,
$\mu_{21}=\\{-0.4,0.2,0.0,0.2\\},\mu_{22}=\\{-0.4,0.2,0.0,0.2\\},\sigma_{1}=\\{0.4,0.6,0.8,1.0,1.2\\},\sigma_{21}=\sigma_{22}=\\{0.05,0.15,0.25,0.35,0.45\\}$
and $\sigma_{\eta}=\\{0.02,0.04,0.10,0.20,1\\}.$

As in the single-agent setting, the standardized doses
$\tilde{d}_{j},\tilde{s}_{l}$ were calibrated in terms of the difference in
the probability of DLE between the neighbouring doses. Specifically,
$\hat{p}_{j}^{(0)}(d_{j})=\hat{p}_{0}^{(0)}+\nu_{d}\times j,$ and
$\hat{p}_{j}^{(0)}(s_{l})=\hat{p}_{0}^{(0)}+\nu_{s}\times j$, $l,j=1,2,3$
where $\nu_{d},\nu_{s}$ are the differences in the toxicity probabilities
between doses of the first and second agents, respectively. The following
values of differences were tried $\nu_{d}=\\{0.025,0.05,0.075,0.10,0.125\\}$,
$\nu_{s}=\\{0.025,0.05,0.075,0.10,0.125\\}$. Finally, we fix $c_{\rm
overdose}=0.25$ for the overdose control constraint.

Using 500 simulations under each scenario and each combination of
hyperparameters, the values
$\mu_{21}=0.0,\mu_{22}=0.0,\sigma_{1}=0.6,\sigma_{21}=\sigma_{22}=0.25$,
$\sigma_{\eta}=0.10,\nu_{d}=\nu_{s}=0.075$ were found to maximize the
geometric mean of the proportion of correct combination selections across the
scenarios.

###  5.3 Results

For the 16 scenarios considered in the simulation study, the percentage of
simulations recommending any desirable dose combination, the percent of
simulations recommending all correct dose combinations and the mean total
sample size are presented in Figure 4, with further detail on individual dose
combination recommendations given in Table 5. The type I error rate, that is
the percentage of simulations recommending any dose combination in Scenario
0-0, is 12.1%. In the extension from monotherapies to dual agent therapies,
some similar patterns are maintained in the results although there are some
notable differences.

Figure 4: Percentage of 10,000 simulations that recommend all desirable dose
combinations (left), the percentage of simulations that recommend any
desirable dose combination (centre) and average total sample size. Note that
only 7 out of 16 efficacy/safety scenarios contain a desirable dose
combination.

It can be seen in Figure 4 that the spread of powers across scenarios for dose
combinations is larger than for monotherapies, for both the selection of all
and any desirable dose combinations. The maximum power of 81.7% to select any
desirable dose combination is achieved in Scenario 0-2, where all doses are
safe and there is a steep monotonic relationship within agent for efficacy.
Even though there is an extra desirable dose combination in Scenario 0-3, we
observe a slightly reduced power since the most efficacious dose combination
has a lower hazard ratio. The power to recommend all desirable doses ranges up
to 61% in Scenario 0-1 where only one dose is desirable, with the lowest in
Scenarios 0-3 and 2-3 where the desirable dose combinations are across
separate agents’ dose escalation (i.e. to select all desirable doses requires
a de-escalation in one agent and then escalation in the other).

Across scenarios, the mean total sample size is 75, ranging from 42-89, a
narrower range than for the single agent. However similar to the single agent,
the smaller sample sizes correspond to scenarios where all dose combinations
are unsafe and therefore the trial stopped early for safety. When this is not
the case, there is little variation across scenarios in terms of mean total
sample size.

Table 5 shows further details of dose recommendations in the simulations.
Especially of note is the emphasis on recommendations of acceptable doses. For
example in Scenario 2-3 where the power to detect desirable dose combinations
is low, a large proportion of simulations also recommend an acceptable dose
combination. It is also clear that inefficacious and/or unsafe doses are
rarely recommended.

Table 5: Percentage of 10,000 simulations where each dose combination is
recommended. Desirable dose combinations are highlighted in bold and
Acceptable dose combinations are highlighted in italics. Note that these may
sum to more than 100% for each scenario as more than one dose combination can
be recommended simultaneously.

Efficacy Scenario 0 1 2 3 0 0.2 0.8 0.2 2.5 0.1 5.1 0.2 5.0 0.8 2.9 2.7 17.3
2.5 25.2 5.0 25.1 0.3 7.9 1.6 61.4 1.4 74.4 2.2 61.3 1 1.1 3.0 1.2 9.5 1.2
17.3 1.3 16.6 5.5 5.3 15.2 23.3 15.1 31.0 28.1 30.7 Safety 0.5 0.5 1.0 1.7 1.2
1.9 1.5 1.7 Scenario 2 2.8 3.2 3.2 9.3 3.1 15.8 3.1 15.8 5.7 4.0 14.7 14.8
15.0 19.3 26.2 19.0 2.1 0.7 8.6 3.6 8.2 4.1 11.1 3.6 3 6.7 1.1 6.5 2.2 6.6 2.8
6.6 2.9 1.0 0.1 2.1 0.2 2.1 0.2 2.9 0.2 0.0 0.0 0.1 0.0 0.1 0.0 0.1 0.0

##  6 Discussion

We introduce and evaluate the statistical design of the AGILE platform which
seeks to quickly establish safe doses and potential for efficacy. We find that
the design can identify potential treatments with good accuracy and show that
the approach is easily extended to combinations of treatments.

The design uses a recently proposed randomized dose-finding design to ensure
that differences between symptoms of COVID can be distinguished from side-
effects of the investigated treatment while a very simple Bayesian model is
used to capture the potential efficacy of the treatments. The latter is in
line with the objective of the trial: make reliable decisions about potential
quickly, rather than using more complex methods that allow more precise
estimation. At the same time this approach guaranteed that the whole platform
structure could be simulated quickly to enable the study design to be fixed
quickly.

The design has been constructed in a flexible manner using a time-to-event
outcome and we based our simulations on time-to-improvement $-$ an endpoint
that has been shown recently to be a highly powered and relatively easy to
collect endpoint (Dodd et al., 2020). The platform has, however, been
constructed to also be able to investigate mild disease in which case a
primary endpoint used would be time-to-negative viral titres in nose and/or
throat swab. Provided that the event rate in this setting is the same, we
expect that the performance reported here will be similar.

## Acknowledgement

T Jaki received funding from UK Medical Research Council (MC_UU_0002/14).
“This report is independent research supported by the National Institute for
Health Research (NIHR Advanced Fellowship, Dr Pavel Mozgunov, NIHR300576; and
Prof Jaki’s Senior Research Fellowship, NIHR-SRF-2015-08-001). The views
expressed in this publication are those of the authors and not necessarily
those of the NHS, the National Institute for Health Research or the Department
of Health and Social Care (DHCS).

## References

  * Bailey et al. (2009) Bailey, S., Neuenschwander, B., Laird, G. and Branson, M. (2009) A Bayesian case study in oncology phase I combination dose-finding using logistic regression with covariates.  Journal of Biopharmaceutical Statistics, 19, 469–484.
  * Beigel et al. (2020) Beigel, J. H., Tomashek, K. M., Dodd, L. E., Mehta, A. K., Zingman, B. S., Kalil, A. C., Hohmann, E., Chu, H. Y., Luetkemeyer, A., Kline, S. et al. (2020) Remdesivir for the treatment of covid-19—preliminary report.  New England Journal of Medicine.
  * Cao et al. (2020) Cao, B., Wang, Y., Wen, D., Liu, W., Wang, J., Fan, G., Ruan, L., Song, B., Cai, Y., Wei, M. et al. (2020) A trial of lopinavir–ritonavir in adults hospitalized with severe covid-19.  New England Journal of Medicine.
  * Dodd et al. (2020) Dodd, L. E., Follmann, D., Wang, J., Koenig, F., Korn, L. L., Schoergenhofer, C., Proschan, M., Hunsberger, S., Bonnett, T., Makowski, M. et al. (2020) Endpoints for randomized controlled clinical trials for covid-19 treatments.  Clinical Trials, 1740774520939938.
  * Griffiths et al. (2020) Griffiths, G., Fitzgerald, R., Jaki, T., Corkhill, A., Marwood, E., Reynolds, H., Stanton, L., Ewings, S., Condie, S., Wrixon, E. et al. (2020) Agile-accord: A randomized, multicentre, seamless, adaptive phase i/ii platform study to determine the optimal dose, safety and efficacy of multiple candidate agents for the treatment of covid-19: A structured summary of a study protocol for a randomised platform trial.  Trials, 21, 1–3.
  * Hirakawa (2012) Hirakawa, A. (2012) An adaptive dose-finding approach for correlated bivariate binary and continuous outcomes in phase I oncology trials.  Statistics in Medicine, 31, 516–532.
  * Horby et al. (2020a) Horby, P., Lim, W. S., Emberson, J. R., Mafham, M., Bell, J. L., Linsell, L., Staplin, N., Brightling, C., Ustianowski, A., Elmahi, E. et al. (2020a) Dexamethasone in hospitalized patients with covid-19-preliminary report.  The New England Journal of Medicine.
  * Horby et al. (2020b) Horby, P., Mafham, M., Linsell, L., Bell, J. L., Staplin, N., Emberson, J. R., Wiselka, M., Ustianowski, A., Elmahi, E., Prudon, B., Whitehouse, A., Felton, T., Williams, J., Faccenda, J., Underwood, J., Baillie, J. K., Chappell, L., Faust, S. N., Jaki, T., Jeffery, K., Lim, W. S., Montgomery, A., Rowan, K., Tarning, J., Watson, J. A., White, N. J., Juszczak, E., Haynes, R. and Landray, M. J. (2020b) Effect of hydroxychloroquine in hospitalized patients with covid-19: Preliminary results from a multi-centre, randomized, controlled trial.  medRxiv.  URL: https://www.medrxiv.org/content/early/2020/07/15/2020.07.15.20151852.
  * Kalbfleisch (1978) Kalbfleisch, J. D. (1978) Nonparametric Bayesian analysis of survival time data.  Journal of the Royal Statistical Society, Series B, 40, 214–221.
  * Lee and Cheung (2009) Lee, S. M. and Cheung, Y. K. (2009) Model calibration in the continual reassessment method.  Clinical Trials, 6, 227–238.
  * Marshall et al. (2020) Marshall, J. C., Murthy, S., Diaz, J., Adhikari, N., Angus, D. C., Arabi, Y. M., Baillie, K., Bauer, M., Berry, S., Blackwood, B. et al. (2020) A minimal common outcome measure set for covid-19 clinical research.  The Lancet Infectious Diseases.
  * Meyer et al. (2020) Meyer, E. L., Mesenbrink, P., Dunger-Baldauf, C., Fülle, H.-J., Glimm, E., Li, Y., Posch, M. and König, F. (2020) The evolution of master protocol clinical trial designs: A systematic literature review.  Clinical Therapeutics.
  * Mozgunov and Jaki (2019a) Mozgunov, P. and Jaki, T. (2019a) A flexible design for advanced Phase I/II clinical trials with continuous efficacy endpoints.  Biometrical Journal, 61, 1477–1492.
  * Mozgunov and Jaki (2019b) — (2019b) An information theoretic phase i–ii design for molecularly targeted agents that does not require an assumption of monotonicity.  Journal of the Royal Statistical Society: Series C (Applied Statistics), 68, 347–367.
  * Mozgunov et al. (2018) Mozgunov, P., Jaki, T. and Paoletti, X. (2018) A benchmark for dose finding studies with continuous outcomes.  Biostatistics, doi: 10.1093/biostatistics/kxy045.
  * Mozgunov et al. (2019) — (2019) Randomized dose-escalation designs for drug combination cancer trials with immunotherapy.  Journal of Biopharmaceutical Statistics, 29, 359–377.
  * Neuenschwander et al. (2008) Neuenschwander, B., Branson, M. and Gsponer, T. (2008) Critical aspects of the Bayesian approach to phase I cancer trials.  Statistics in Medicine, 27, 2420–2439.
  * Neuenschwander et al. (2015) Neuenschwander, B., Matano, A., Tang, Z., Roychoudhury, S., Wandel, S. and Bailey, S. (2015) A Bayesian industry approach to phase I combination trials in oncology.  In Statistical methods in drug combination studies (eds. W. Zhao and H. Yang), chap. 6, 95–135. Chapman & Hall/CRC Press: Boca Raton, FL.
  * Pallmann et al. (2018) Pallmann, P., Bedding, A. W., Choodari-Oskooei, B., Dimairo, M., Flight, L., Hampson, L. V., Holmes, J., Mander, A. P., Sydes, M. R., Villar, S. S. et al. (2018) Adaptive designs in clinical trials: why use them, and how to run and report them.  BMC Medicine, 16, 1–15.
  * Schmidli et al. (2014) Schmidli, H., Gsteiger, S., Roychoudhury, S., O’Hagan, A., Spiegelhalter, D. and Neuenschwander, B. (2014) Robust meta-analytic-predictive priors in clinical trials with historical control information.  Biometrics, 70, 1023–1032.
  * Sinha et al. (2003) Sinha, D., Ibrahim, J. and Chen, M.-H. (2003) A Bayesian justification of Cox’s partial likelihood.  Biometrika, 90, 629–641.
  * Stallard et al. (2020) Stallard, N., Hampson, L., Benda, N., Brannath, W., Burnett, T., Friede, T., Kimani, P. K., Koenig, F., Krisam, J., Mozgunov, P. et al. (2020) Efficient adaptive designs for clinical trials of interventions for covid-19.  Statistics in Biopharmaceutical Research.
  * Wages and Tait (2015) Wages, N. A. and Tait, C. (2015) Seamless phase I/II adaptive design for oncology trials of molecularly targeted agents.  Journal of Biopharmaceutical Statistics, 25, 903–920.
  * Wang et al. (2020) Wang, Y., Zhang, D., Du, G., Du, R., Zhao, J., Jin, Y., Fu, S., Gao, L., Cheng, Z., Lu, Q. et al. (2020) Remdesivir in adults with severe covid-19: a randomised, double-blind, placebo-controlled, multicentre trial.  The Lancet.


# Cell-Free Massive MIMO for ISAC: Access Point Operation Mode Selection and
Power Control

Mohamed Elfiatoure, Mohammadali Mohammadi, Hien Quoc Ngo, and Michail
Matthaiou
Centre for Wireless Innovation (CWI), Queen’s University Belfast, U.K.
Email:{melfiatoure01, m.mohammadi, hien.ngo<EMAIL_ADDRESS>
###### Abstract

This paper considers a cell-free massive multiple-input multiple-output (MIMO)
integrated sensing and communication (ISAC) system, where distributed MIMO
access points (APs) are used to jointly serve the communication users and
detect the presence of a single target. We investigate the problem of AP
operation mode selection, wherein some APs are dedicated for downlink
communication, while the remaining APs are used for sensing purposes. Closed-
form expressions for the individual spectral efficiency (SE) and mainlobe-to-
average-sidelobe ratio (MASR) are derived, which are respectively utilized to
assess the communication and sensing performances. Accordingly, a max–min
fairness problem is formulated and solved, where the minimum SE of the users
is maximized, subject to the per-AP power constraints as well as sensing MASR
constraint. Our numerical results show that the proposed AP operation mode
selection with power control can significantly improve the communication
performance for given sensing requirements.

††footnotetext: This project has received funding from the European Research
Council (ERC) under the European Union’s Horizon 2020 research and innovation
programme (grant agreement No. 101001331). The work of H. Q. Ngo was supported
by the U.K. Research and Innovation Future Leaders Fellowships under Grant
MR/X010635/1.

##  I Introduction

ISAC has recently been envisioned as a key enabling technology for future
wireless networks, aiming to efficiently utilize the congested resources for
both communication and sensing [1, 2]. The radar bands set aside for sensing
can be harnessed for wireless communication operation, enabling the
implementation of high data-rate applications. To unify the radar and
communication operations, two well-known designs, namely separated and co-
located systems, were introduced in [3, 4, 5] and [6, 7], respectively. The
former utilizes different devices, operating over the same frequency band, for
communication and sensing, while in the latter a single device acts as radar
and communication base station (BS) by simultaneously communicating with
multiple downlink users and detecting radar targets.

The main driving force behind the transition from the separated design to a
co-located design was to reduce the complexity induced by side-information
exchange among the radar and communication devices [7]. However, co-located
design with a MIMO BS often suffers from a fairness problem, since the cell-
boundary users are subject to inter-cell interference and significant power
decay over long distances. The key feature of massive MIMO technology, i.e.,
inter/intra-cell interference suppression, revitalizes the interest towards
separated design with multiple communication and radar devices to implement
distributed ISAC architectures. In this context, cell-free massive MIMO with
distributed MIMO APs can be exploited to support ISAC. In cell-free massive
MIMO, all users are coherently served by all APs over the same time-frequency
band. Each AP is connected to a central processing unit (CPU) via fronthaul
links, and the CPU is responsible for coordination [8, 9].

The integration of ISAC into cell-free massive MIMO networks, has been
recently investigated in [10, 11]. Specifically, Behdad et al. [10] studied a
cell-free massive MIMO ISAC system, consisting of a fixed number of transmit
and receive APs. Users are served by the transmit APs and, at the same time,
the transmitted signals are used for sensing to detect the presence of a
target in a certain location. The reflected signals are received at the
receive APs and then processed at the CPU. The authors proposed a power
allocation algorithm to maximize the sensing signal-to-noise ratio (SNR) under
signal-to-interference-plus-noise ratio (SINR) constraints at the user.
Demirhan et al. [11] studied the sensing and communication beamforming design
problem in cell-free massive MIMO ISAC systems, where a joint beamforming
design was proposed to maximize the sensing SNR, while satisfying the
communication SINR constraints.

Different from the above-mentioned works [10, 11], where the AP operation
modes are fixed, we consider a novel cell-free massive MIMO ISAC network with
dynamic AP operation mode selection. The APs’ operation mode is designed to
maximize the minimum SE of the downlink users, while satisfying the sensing
requirement to detect the presence of a single target in a certain location.
Relying on the long-term channel state information (CSI), the APs are divided
into communication APs (C-APs) and sensing APs (S-AP) to support downlink
communication and sensing operations simultaneously. The main contributions of
our paper can be summarized as follows:

  * •

By leveraging the use-and-then-forget strategy, we derive closed-form
expressions for the downlink SE and MRSR to evaluate the performance of the
communication and sensing operation, respectively. Then, we formulate the
problem of joint AP operation mode selection and power control, considering
per-AP power constraints and a MASR constraint for target detection.

  * •

We propose a greedy algorithm for AP operation mode selection. Accordingly, an
alternating optimization (AO) algorithm is developed to handle the coupling
between the C-AP and R-AP power control coefficients’ design.

  * •

Numerical results show that our proposed greedy AP operation mode selection
with optimal power control (GAP-OPC) significantly improves the SE performance
of the downlink users for given MASR, compared to the greedy/random operation
mode selection with no power control (GAP/RAP-NPC) benchmarks.

Notation: We use bold lower case letters to denote vectors. The superscripts
$(\cdot)^{*}$ and $(\cdot)^{T}$ stand for the conjugate and transpose,
respectively; $\mathbf{I}_{N}$ denotes the $N\times N$ identity matrix. A zero
mean circular symmetric complex Gaussian distribution having variance
$\sigma^{2}$ is denoted by $\mathcal{CN}(0,\sigma^{2})$. Finally,
$\mathbb{E}\\{\cdot\\}$ denotes the statistical expectation.

##  II System Model

We consider a cell-free massive MIMO ISAC system under time division duplex
operation, where $M$ APs serve $K_{d}$ users in the downlink, while radiating
probing signals to a target direction for radar sensing. Each user is equipped
with one single antenna, while each AP is equipped with $N$ antennas. All APs
and users operate as half-duplex devices. For notational simplicity, we define
the sets $\mathcal{M}\triangleq\\{1,\ldots,M\\}$ and
$\mathcal{K}_{d}\triangleq\\{1,\dots,K_{d}\\}$ as the collections of indices
of the APs and users, respectively. As shown in Fig. 1, downlink communication
as well as target detection take place simultaneously and over the same
frequency band. The AP operation mode selection approach is designed according
to the network requirements, determining whether an AP is dedicated to
information transmission or radar sensing. The users receive information from
a group of the APs, termed as C-APs, while the remaining APs, termed as S-APs,
are used for target detection.

###  II-A Channel Model and Uplink Training

We assume a quasi-static channel model, with each channel coherence interval
spanning a duration of $\tau$ symbols. The duration of the training is denoted
as $\tau_{t}$, while the duration of downlink information transfer and target
detection is $(\tau-\tau_{t})$.

For the sensing channel model, we assume there is a line-of-sight (LOS) path
between the target location and each AP, which is a commonly adopted model in
the literature [10, 11]. The LOS channel between AP $m$ and target is given by

$\displaystyle~{}{{\bf g}}_{m}$ $\displaystyle={\bf
a}_{N}(\phi_{m,t}^{a},\phi_{m,t}^{e}),~{}\forall m\in\mathcal{M},$ (1)

where $\phi_{m,t}^{a}$, $\phi_{m,t}^{e}$ denote the azimuth and elevation
angles of departure (AoD) from AP $m$ towards the target. Moreover, the $q$-th
entry of the array response vector ${\bf
a}_{N}(\phi_{m,t}^{a},\phi_{m,t}^{e})\in\mathbb{C}^{N\times 1}$, is given by

$\displaystyle~{}[{\bf a}_{N}(\phi_{m,t}^{a},\phi_{m,t}^{e})]_{q}$
$\displaystyle\\!=\\!\frac{1}{\sqrt{N}}\exp\Big{(}j2\pi\frac{d}{\lambda}(q\\!-\\!1)\sin\phi_{m,t}^{e}\sin\phi_{m,t}^{a}\Big{)},$
(2)

where $d$ and $\lambda$ denote the AP antenna spacing and carrier wavelength,
respectively.

The channel vector between the $m$-th AP and $k$-th user is modeled as
$\mathbf{g}_{mk}=\sqrt{\beta_{mk}}{\bf h}_{mk}$, where $\beta_{mk}$ is the
large scale fading coefficients, and ${\bf h}_{mk}\in\mathbb{C}^{N\times 1}$
is the small-scale fading vector, whose elements are independent and
identically distributed $\mathcal{CN}(0,1)$ random variables [9].

An uplink training process is implemented to acquire the local CSI between
each AP and all users. In each coherence block of length $\tau$, all users are
assumed to transmit their pairwisely orthogonal pilot sequence of length
$\tau_{t}$ to all APs, which requires $\tau_{t}\geq K_{d}$. At AP $m$, ${\bf
g}_{mk}$ is estimated by using the received pilot signals and applying the
minimum mean-square error (MMSE) estimation technique. By following [9], the
MMSE estimate $\hat{{\bf g}}_{mk}$ of ${\bf g}_{mk}$ is obtained as $\hat{{\bf
g}}_{mk}\sim\mathcal{CN}\left(\boldsymbol{0},\gamma_{mk}{\bf I}_{N}\right)$,
where

$\displaystyle\gamma_{mk}$
$\displaystyle=\frac{\tau_{t}\rho_{t}\beta_{mk}^{2}}{\tau_{t}\rho_{t}\beta_{mk}+1},$
(3)

while $\rho_{t}$ represents the normalized transmit power of each pilot
symbol.

Figure 1: Cell-free massive MIMO ISAC system.

###  II-B Data and Probing Signal Transmission

AP operation mode selection is performed by considering large-scale fading
effects and relying on the statistical CSI, obtained during the training
phase. The binary variables used to indicate the operation mode for each AP
$m$ are defined as

$\displaystyle a_{m}=\begin{cases}1,&\text{if AP $m$ operates as C-AP}\\\
0,&\text{if AP $m$ operates as S-AP.}\end{cases}$ (4)

The transmission phase comprises information transmission from C-APs to users
and probing signal transmission from S-APs to the target. Let
$\mathbf{x}_{c,m}$ and $\mathbf{x}_{r,m}$ denote the data and probing signals,
respectively, where $\mathbb{E}\big{\\{}|{\bf x}_{c,m}|^{2}\big{\\}}=1$ and
$\mathbb{E}\big{\\{}|{\bf x}_{r,m}|^{2}\big{\\}}=1$. The signal vector
transmitted from AP $m$ can be expressed as

$\displaystyle\mathbf{x}_{m}={a_{m}}{{\bf x}}_{c,m}+({1-a_{m}}){{\bf
x}_{r,m}}.$ (5)

The power control coefficients at AP $m$ are chosen to satisfy the power
constraint at each S-AP and C-AP, respectively, i.e.,

$\displaystyle{a_{m}}\mathbb{E}\\{\|\mathbf{x}_{c,m}\|^{2}\\}+({1-a_{m}})\mathbb{E}\\{\|\mathbf{x}_{r,m}\|^{2}\\}\leq\rho,$
(6)

where $\rho$ denotes the maximum normalized downlink power.

The transmit signal for communication at the $m$-th C-AP can be expressed as
$\mathbf{x}_{c,m}=\sum_{k\in\mathcal{K}_{d}}\sqrt{{\eta}_{mk}\rho}{\bf
t}_{mk}^{\mathrm{Com}}x_{c,k}$, where $\eta_{mk}$ represents the downlink
power control coefficient at the $m$-th C-AP during communication, while ${\bf
t}_{mk}^{\mathrm{Com}}\in\mathbb{C}^{N\times 1}$ and $x_{c,k}$ denote the
precoding vector and intended signal for user $k$, respectively.

Moreover, the probing signal transmitted by the $m$-th S-AP can be expressed
as ${\bf x}_{r,m}=\sqrt{\eta_{m}\rho}{\bf t}_{m}^{\mathrm{Sen}}x_{r}$, where
$\eta_{m}$ denotes the power control coefficient at S-AP $m$, while ${\bf
t}_{m}^{\mathrm{Sen}}\in\mathbb{C}^{N\times 1}$ denotes the beamforming vector
for sensing and $x_{r}$ is the sensing symbol.

We note that in the absence of communication users, conjugate sensing
beamforming solution becomes optimal, as it directly maximizes the sensing SNR
[10, 11]. Therefore, to detect the presence of a target in a certain location,
we design the sensing beamforming vector at S-AP $m$ as

$\displaystyle~{}{\bf t}_{m}^{\mathrm{Sen}}={\bf
a}_{N}(\phi_{m,t}^{a},\phi_{m,t}^{e}).$ (7)

Furthermore, the conjugate precoder is employed at the C-AP, particularly due
to its advantages including low computational complexity, ease of analysis,
and reasonable performance, as shown in [12, 13]. Hence, ${\bf
t}_{mk}^{\mathrm{Com}}$ is given by

$\displaystyle~{}{\bf t}_{mk}^{\mathrm{Com}}=\hat{{\bf g}}_{mk}^{*}.$ (8)

###  II-C Sensing Operation and MASR

For a given channel realization, the average spatial power pattern for sensing
is defined as

$\displaystyle~{}P^{\mathrm{ave}}(\phi_{m,t}^{a},\phi_{m,t}^{e})$
$\displaystyle=\mathbb{E}\Bigg{\\{}\sum_{m\in\mathcal{M}}|{\bf
a}_{N}^{H}(\phi_{m,t}^{a},\phi_{m,t}^{e}){\bf x}_{m}|^{2}\Bigg{\\}}$
$\displaystyle=\rho\sum_{m\in\mathcal{M}}{a_{m}}\mathbb{E}\Bigg{\\{}\bigg{|}\sum_{k\in\mathcal{K}_{d}}\sqrt{{\eta}_{mk}}{\bf
a}_{N}^{H}(\phi_{m,t}^{a},\phi_{m,t}^{e}){\bf
t}_{mk}^{\mathrm{Com}}\bigg{|}^{2}\Bigg{\\}}$
$\displaystyle+\rho\sum_{m\in\mathcal{M}}(1-a_{m})\eta_{m}\mathbb{E}\Big{\\{}\big{|}{\bf
a}_{N}^{H}(\phi_{m,t}^{a},\phi_{m,t}^{e}){\bf
t}_{m}^{\mathrm{Sen}}\big{|}^{2}\Big{\\}},$ (9)

where the expectation is taken over the transmitted signals, assuming that the
information signal and probing signals are independent zero-mean Gaussian
distributed.

###### Proposition 1

The average spatial power pattern for sensing is given by
$P^{\mathrm{ave}}(\phi_{m,t}^{a},\phi_{m,t}^{e})=P^{\mathrm{ave}}_{\mathrm{Com}}(\phi_{m,t}^{a},\phi_{m,t}^{e})+P^{\mathrm{ave}}_{\mathrm{Sen}}(\phi_{m,t}^{a},\phi_{m,t}^{e})$,
where

$\displaystyle P^{\mathrm{ave}}_{\mathrm{Com}}(\phi_{m,t}^{a},\phi_{m,t}^{e})$
$\displaystyle=\rho\sum_{m\in\mathcal{M}}\sum_{k\in\mathcal{K}_{d}}a_{m}{\eta}_{mk}\gamma_{mk},$
(10a) $\displaystyle
P^{\mathrm{ave}}_{\mathrm{Sen}}(\phi_{m,t}^{a},\phi_{m,t}^{e})$
$\displaystyle=\rho\sum_{m\in\mathcal{M}}(1-a_{m})\eta_{m}.$ (10b)

Proof: By invoking (7) and (8), and then taking the expectation of (II-C) over
${\bf t}_{mk}^{\mathrm{Com}}$, the desired result is obtained. $\hskip
40.00006pt\blacksquare$

We would like $P^{\mathrm{ave}}{\mathrm{Com}}(\phi{m,t}^{a},\phi_{m,t}^{e})$,
$\forall\phi_{m,t}^{a},\phi_{m,t}^{e}$, to be as small as possible to confine
the pattern distortion. For illuminating a target angle
$(\phi_{m,t}^{a},\phi_{m,t}^{e})$, it is desirable that the mainlobe level
$P^{\mathrm{ave}}_{\mathrm{Sen}}(\phi_{m,t}^{a},\phi_{m,t}^{e})$ is higher
than $P^{\mathrm{ave}}_{\mathrm{Com}}(\phi_{m,t}^{a},\phi_{m,t}^{e})$ by a
certain minimum sensing level $\kappa$, which is referred to as the MASR:

$\displaystyle~{}\mathrm{MASR}({\bf
a},\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}})$
$\displaystyle=\frac{P^{\mathrm{ave}}_{\mathrm{Sen}}(\phi_{m,t}^{a},\phi_{m,t}^{e})}{P^{\mathrm{ave}}_{\mathrm{Com}}(\phi_{m,t}^{a},\phi_{m,t}^{e})}$
(11)
$\displaystyle=\frac{\sum_{m\in\mathcal{M}}(1-a_{m})\eta_{m}}{\sum_{m\in\mathcal{M}}\sum_{k\in\mathcal{K}_{d}}a_{m}{\eta}_{mk}\gamma_{mk}}\geq\kappa,$

where ${\bf a}=\\{a_{1},\ldots,a_{M}\\}$,
$\boldsymbol{\eta}^{\mathrm{Com}}=\\{\eta_{m1},\ldots,\eta_{mK_{d}}\\}$,
$\forall m\in\mathcal{M}$, and
$\boldsymbol{\eta}^{\mathrm{Sen}}=\\{\eta_{1},\ldots,\eta_{M}\\}$.

###  II-D Communication Operation and Downlink SE

The received signal at $k$-th user can be expressed as

$\displaystyle{y}_{k}$
$\displaystyle={\sum_{m\in\mathcal{M}}a_{m}\sqrt{\rho\eta_{mk}}\mathbf{g}^{T}_{mk}{\bf
t}_{mk}^{\mathrm{Com}}x_{c,k}}$ $\displaystyle\hskip
20.00003pt+{\sum_{m\in\mathcal{M}}\sum_{k^{\prime}\in\mathcal{K}_{d}\setminus
k}a_{m}\sqrt{\rho\eta_{mk^{\prime}}}\mathbf{g}^{T}_{mk}{\bf
t}_{mk^{\prime}}^{\mathrm{Com}}x_{c,k^{\prime}}}$ $\displaystyle\hskip
20.00003pt+{\sum_{m\in\mathcal{M}}\left(1-a_{m}\right)\sqrt{\rho\eta_{m}}{\bf
g}^{T}_{mk}{\bf t}_{m}^{\mathrm{Sen}}x_{r}}+{n}_{k},$ (12)

where the second term is the inter-user interference, the third term
represents the interference from S-APs, and
$n_{k}\sim\mathcal{CN}(0,\sigma^{2}_{n})$ denotes the additive white Gaussian
noise at the user $k$.

###### Proposition 2

With conjugate precoding at the APs for downlink communication, the achievable
downlink SE of user $k$, can be expressed as
$\mathrm{SE}_{k}=\Big{(}1-\frac{\tau_{\mathrm{p}}}{\tau}\Big{)}\log_{2}\big{(}1+\mathrm{SINR}_{k}({\bf
a},\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}})\big{)}$,
where $\mathrm{SINR}_{k}({\bf
a},\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}})$ is
given by (13) at the top of the next page.

$\displaystyle\mathrm{SINR}_{k}({\bf
a},\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}})=\\!\frac{\rho
N^{2}\big{(}\sum_{m\in\mathcal{M}}\\!a_{m}\eta_{mk}^{1/2}\gamma_{mk}\big{)}^{2}}{\rho
N\\!\\!\sum_{m\in\mathcal{M}}\\!\sum_{k^{\prime}\in\mathcal{K}_{d}}\\!a_{m}\eta_{mk^{\prime}}\gamma_{mk^{\prime}}\beta_{mk}\\!+\\!\rho\\!\sum_{m\in\mathcal{M}}(1\\!-\\!a_{m}){\eta_{m}}\beta_{mk}\\!+\\!1}\\!.$
(13)

Proof: See Appendix A. $\hskip 140.00021pt\blacksquare$

##  III Proposed Design Problem and Solution

In this section, we formulate and solve the AP mode selection to maximize the
minimum SE. More specifically, we aim to optimize the AP operation mode
selection vector (${\bf a}$) and power control coefficients
($\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}}$) to
maximize the minimum per-user SE subject to a prescribed MASR level for the
target detection and transmit power constraints at the APs. The optimization
problem is then formulated as

$\displaystyle\text{{(P1):}}~{}\underset{{\bf
a},\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}}}{\max}$
$\displaystyle\min_{k\in\mathcal{K}_{d}}\,\,\hskip 5.0pt\mathrm{SINR}_{k}({\bf
a},\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}})$ (14a)
$\displaystyle\mathrm{s.t.}\,\,\hskip 5.0pt$ $\displaystyle\mathrm{MASR}({\bf
a},\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}})\geq\kappa,$
(14b) $\displaystyle
a_{m}\sum_{k\in\mathcal{K}_{d}}{\eta_{mk}}\gamma_{mk}\leq\frac{1}{N},~{}\forall
m\in\mathcal{M},$ (14c) $\displaystyle{\eta_{m}}\leq 1-a_{m},~{}\forall
m\in\mathcal{M},$ (14d) $\displaystyle a_{m}\in\\{0,1\\}.$ (14e)

Problem (P1) is a challenging combinatorial problem. Therefore, for AP
operation mode selection, we only focus on a heuristic greedy method which
simplifies the computation, while providing a significantly successful
monitoring performance gain.

###  III-A AP Operation Mode Selection with Fixed Power Control

Let $\mathcal{A}_{\mathrm{Sen}}$ and $\mathcal{A}_{\mathrm{Com}}$ denote the
sets containing the indices of APs operating as radar, i.e., APs with
$a_{m}=0$, and APs operating in communication mode, i.e., APs with $a_{m}=1$,
respectively. In addition,
$\mathrm{MASR}(\mathcal{A}_{\mathrm{Sen}},\mathcal{A}_{\mathrm{Com}})$ and
$\mathrm{SINR}_{k}(\mathcal{A}_{\mathrm{Sen}},\mathcal{A}_{\mathrm{Com}})$
underline the dependence of the sensing MASR and received SINR of the $k$-th
user on the different choices of AP mode selections. Our greedy algorithm of
AP mode selection is shown in Algorithm 1. To guarantee the sensing MASR
requirement, all APs are initially assigned for sensing operation, i.e.,
$\mathcal{A}_{\mathrm{Sen}}=\mathcal{M}$ and
$\mathcal{A}_{\mathrm{Com}}=\emptyset$. Then, in each iteration, one AP
switches into communication operation mode for maximizing the minimum of SE
(or equivalently SINR in (13)), while the minimum MARS required for target
sensing is guaranteed. This process continues until there is no more
improvement in the minimum SINR among all users.

Algorithm 1 Greedy AP Operation Mode Selection

1: Initialize: Set $\mathcal{A}_{\mathrm{Com}}=\emptyset$ and
$\mathcal{A}_{\mathrm{Sen}}=\mathcal{M}$. Set iteration index $i=0$.

2: Calculate
$\Pi^{\star}[i]=\min_{k\in\mathcal{K}_{d}}\mathrm{SE}_{k}({A}_{\mathrm{Sen}},\mathcal{A}_{\mathrm{Com}})$

3: repeat

4: for all $m\in\mathcal{A}_{\mathrm{Sen}}$ do

5: Set $\mathcal{A}_{s}=\mathcal{A}_{\mathrm{Sen}}\setminus m$.

6: if $\mathrm{MASR}(\mathcal{A}_{s},\mathcal{A}_{\mathrm{Com}}\bigcup
m)\\!\geq\\!\kappa$ then

7: Calculate
$\Pi_{m}=\min_{k\in\mathcal{K}_{d}}\mathrm{SINR}_{k}(\mathcal{A}_{s},\mathcal{A}_{\mathrm{Com}}\bigcup
m)$

8: else

9: Set $\Pi_{m}=0$

10: end if

11: end for

12: Set
$\Pi^{\star}[i+1]=\underset{m\in\mathcal{A}_{\mathrm{Com}}}{\max}\,\,\Pi_{m}$

13: $e=|\Pi^{\star}[i+1]-\Pi^{\star}[i]|$

14: if $e\geq e_{\min}$ then

15: Update
$\mathcal{A}_{\mathrm{Com}}\\!=\\!\\{\mathcal{A}_{\mathrm{Com}}\bigcup
m^{\star}\\}$ and
$\mathcal{A}_{\mathrm{Sen}}\\!=\\!\mathcal{A}_{\mathrm{Sen}}\\!\setminus\\!m^{\star}$

16: end if

17: Set $i=i+1$

18: until $e<e_{\min}$

19: return $\mathcal{A}_{\mathrm{Sen}}$ and $\mathcal{A}_{\mathrm{Com}}$,
i.e., the indices of APs operating in radar mode and communication mode,
respectively.

###  III-B Power Control

For a given AP mode selection, the optimization problem (14) reduces to the
power control problem, given by

$\displaystyle\text{{(P2):}}~{}\underset{\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}}}{\max}\,\,\hskip
5.0pt$ $\displaystyle\min_{k\in\mathcal{K}_{d}}\,\,\hskip
5.0pt\mathrm{SINR}_{k}(\boldsymbol{\eta}^{\mathrm{Com}},\boldsymbol{\eta}^{\mathrm{Sen}})$
(15a) $\displaystyle\mathrm{s.t.}\,\,\hskip
5.0pt\eqref{eq:MASR:ct1}-\eqref{eq:am:binary}.$ (15b)

Problem (P2) is a non-convex optimization problem due the to non-convex
objective function and constraints. Since the variables
$\boldsymbol{\eta}^{\mathrm{Com}}$ and $\boldsymbol{\eta}^{\mathrm{Sen}}$ are
coupled in both the objective and MASR constraint, it is difficult to
simultaneously optimize them. Therefore, we propose an AO algorithm to jointly
optimize $\boldsymbol{\eta}^{\mathrm{Com}}$ and
$\boldsymbol{\eta}^{\mathrm{Sen}}$ in two sub-problems.

Firstly, for a given $\boldsymbol{\eta}^{\mathrm{Sen}}$, we formulate the sub-
problem for optimizing $\boldsymbol{\eta}^{\mathrm{Com}}$ as

$\displaystyle\text{{(P2-1):}}~{}\underset{\boldsymbol{\eta}^{\mathrm{Com}}}{\max}\,\,\hskip
5.0pt$ $\displaystyle\min_{k\in\mathcal{K}_{d}}\,\,\hskip
5.0pt\mathrm{SINR}_{k}(\boldsymbol{\eta}^{\mathrm{Com}})$ (16a)
$\displaystyle\mathrm{s.t.}\,\,\hskip
5.0pt\eqref{eq:MASR:ct1},\eqref{eq:etamk:ct1}.$ (16b)

By introducing the slack variables $\theta_{mk}=\eta_{mk}^{\frac{1}{2}}$ and
$\upsilon_{m}$, we reformulate (P2-1) as

$\displaystyle\text{{(P2-2):}}~{}\underset{\boldsymbol{\eta}^{\mathrm{Com}},t}{\max}\,\,\hskip
5.0ptt$ (17a) $\displaystyle\mathrm{s.t.}\,\,\hskip
5.0pt\\!\frac{\big{(}\sum_{m\in\mathcal{M}}\\!a_{m}\theta_{mk}\gamma_{mk}\big{)}^{2}}{\frac{1}{N}\\!\\!\sum_{m\in\mathcal{M}}a_{m}\beta_{mk}\upsilon_{m}^{2}\\!+\varphi_{k}}\geq
t,~{}\forall k\in\mathcal{K}_{d}\\!$ (17b)
$\displaystyle\frac{\sum_{m\in\mathcal{M}}(1-a_{m})\eta_{m}}{\sum_{m\in\mathcal{M}}a_{m}\upsilon_{m}^{2}}\geq\kappa,$
(17c)
$\displaystyle\sum_{k^{\prime}\in\mathcal{K}_{d}}\\!a_{m}\theta_{mk^{\prime}}^{2}\leq\upsilon_{m}^{2},~{}\forall
m\in\mathcal{M}$ (17d) $\displaystyle 0\leq
a_{m}\upsilon_{m}^{2}\leq\frac{1}{N},~{}\forall m\in\mathcal{M},$ (17e)
$\displaystyle\theta_{mk}\geq 0,\forall m\in\mathcal{M},\forall
k\in\mathcal{K}_{d},$ (17f)

where
$\varphi_{k}\overset{\Delta}{=}\frac{1}{N^{2}}\\!\sum_{m\in\mathcal{M}}(1\\!-\\!a_{m}){\eta_{m}}\beta_{mk}\\!+\frac{1}{\rho
N^{2}}$ and and
$\upsilon_{m}^{2}\overset{\Delta}{=}\sum_{k\in\mathcal{K}_{d}}\\!a_{m}\eta_{mk}$.
The equivalence between (17) and (16) follows directly from the fact that the
second constraint in (17) holds with equality at the optimum. Problem (P2-2)
can be reformulated as a second-order cone program (SOCP). More precisely, for
given $t$, we have

$\displaystyle\text{{(P2-3):}}~{}\underset{\boldsymbol{\eta}^{\mathrm{Com}}}{\max}\,\,\hskip
5.0ptt$ (18) $\displaystyle\mathrm{s.t.}\,\,\hskip 5.0pt\|{\bf
v}_{k}\|\leq\frac{1}{\sqrt{t}}\sum_{m\in\mathcal{M}}\\!a_{m}\theta_{mk}\gamma_{mk},~{}\forall
k\in\mathcal{K}_{d},$ (19)
$\displaystyle~{}\eqref{eq:ct2:P2-2}-~{}\eqref{eq:ct5:P2-2},$ (20)

where ${\bf v}_{k}=[{\bf v}_{k1}^{T},\sqrt{\varphi_{k}}]^{T}$, with ${\bf
v}_{k1}=\big{[}\sqrt{\frac{\beta_{1k}}{N}}\upsilon_{1},\ldots,\sqrt{\frac{\beta_{Mk}}{N}}\upsilon_{M}\big{]}^{T}$.
The first constraint represents a second order cone and thus (P2-3) is a
standard SOCP, which is a convex problem. The bisection search method is
exploited to find the optimal solution, in each step solving a sequence of
convex feasibility problem. This bisection based search method is summarized
in Algorithm 2.

Algorithm 2 Bisection Method for Power Control

1: Initialization of $t_{\min}$ and $t_{\max}$, where $t_{\min}$ and
$t_{\max}$ define a range of relevant values of the objective function in
(16). Initial line-search accuracy $\epsilon$.

2: repeat

3: Set $t:=\frac{t_{\min}+t_{\max}}{2}$. Solve the following convex
feasibility program  $\displaystyle\begin{cases}&\|{\bf
v}_{k}\|\leq\frac{1}{\sqrt{t}}\big{(}\sum_{m\in\mathcal{M}}\\!a_{m}\theta_{mk}\gamma_{mk}\big{)},~{}\forall
k\in\mathcal{K}_{d}\\!.\\\ &\hskip
5.0pt{\sum_{m\in\mathcal{M}}a_{m}\upsilon_{m}^{2}}\leq\frac{1}{\kappa}{\sum_{m\in\mathcal{M}}(1-a_{m})\eta_{m}},\\\
&\sum_{k^{\prime}\in\mathcal{K}_{d}}\\!a_{m}\theta_{mk^{\prime}}^{2}\leq\upsilon_{m}^{2},~{}\forall
m\in\mathcal{M}\\\ &0\leq a_{m}\upsilon_{m}\leq\frac{1}{\sqrt{N}},~{}\forall
m\in\mathcal{M},\\\ &\theta_{mk}\geq 0,\forall m\in\mathcal{M},\forall
k\in\mathcal{K}_{d},\end{cases}$ (21)

4: If problem (21) is feasible, then set $t_{\min}:=t$, else set
$t_{\max}:=t$.

5: until $t_{\max}-t_{\min}<\epsilon$.

Secondly, when $\boldsymbol{\eta}^{\mathrm{Com}}$ is fixed, the sub-problem
for optimizing $\boldsymbol{\eta}^{\mathrm{Sen}}$ can be expressed as

$\displaystyle\text{{(P2-4):}}~{}\underset{\boldsymbol{\eta}^{\mathrm{Sen}}}{\max}\,\,\hskip
5.0pt$ $\displaystyle\min_{k\in\mathcal{K}_{d}}\,\,\hskip
5.0pt\mathrm{SINR}_{k}(\boldsymbol{\eta}^{\mathrm{Sen}})$ (22a)
$\displaystyle\mathrm{s.t.}\,\,\hskip
5.0pt\eqref{eq:MASR:ct1},\eqref{eq:etam:ct1}.$ (22b)

By introducing a new slack variable $\varrho$, we can reformulate the
optimization problem as

$\displaystyle\text{{(P2-5):}}~{}\underset{\boldsymbol{\eta}^{\mathrm{Sen}},\varrho}{\max}\,\,\hskip
5.0pt\varrho$ (23a) $\displaystyle\mathrm{s.t.}\,\,\hskip
5.0pt\\!\frac{\\!\sum_{m\in\mathcal{M}}(1\\!-\\!a_{m}){\eta_{m}}\beta_{mk}\\!+\phi_{k}}{\big{(}N\sum_{m\in\mathcal{M}}\\!a_{m}\theta_{mk}\gamma_{mk}\big{)}^{2}}\leq\frac{1}{\varrho},~{}\forall
k\in\mathcal{K}_{d}\\!$ (23b)
$\displaystyle\sum_{m\in\mathcal{M}}(1-a_{m})\eta_{m}\geq\kappa\sum_{m\in\mathcal{M}}a_{m}\upsilon_{m}^{2},$
(23c) $\displaystyle{\eta_{m}}\leq 1-a_{m},~{}\forall m\in\mathcal{M},$ (23d)

where
$\phi_{k}\overset{\Delta}{=}\frac{1}{N}\\!\\!\sum_{m\in\mathcal{M}}a_{m}\beta_{mk}\upsilon_{m}^{2}\\!+\frac{1}{\rho
N^{2}}$. Now, for a fixed $\varrho$, all inequalities involved in (P2-5) are
linear, hence the solution to the optimization problem can be obtained by
harnessing a line-search over $\varrho$ to find the maximal feasible value.
Therefore, we can apply the bisection method in Algorithm 2 to solve (23),
where $t_{\min}$, $t_{\max}$, and the feasibility problem (21) are replaced
with $\varrho_{\min}$, $\varrho_{\max}$ and problem (23), respectively. We
summarize the overall AO algorithm in Algorithm 3.

###  III-C Complexity Analysis

Here, we provide the computational complexity of Algorithm 3, which involves a
SOCP problem in (18) and a linear-search problem (23) at each iteration. In
order to solve a SOCP, the iterative bisection search method requires
$\mathcal{O}(n_{v_{1}}^{2}n_{c})$ arithmetic operations, where $n_{v_{1}}$ is
the number of optimization variables and $n_{c}$ is the total number of SOC
constraints [14]. Moreover, the total number of iterations required is
$\log_{2}\big{(}\frac{t_{\max}-t_{\min}}{\epsilon}\big{)}$. In (18), the total
number of variables is $n_{v_{1}}=MK_{d}$ and there are $n_{v}=K_{d}$ SOC
constraints. Therefore, the per-iteration computational complexity for solving
(18) is
$\log_{2}\big{(}\frac{t_{\max}-t_{\min}}{\epsilon}\big{)}\mathcal{O}(MK_{d}^{3})$.
Problem (23) involves $n_{v_{2}}=M$ scalar-value variables and
$n_{c_{2}}=M+K_{d}+1$ linear constraints. According to [15], the per-iteration
cost to solve (23) is
$\mathcal{O}\big{(}(n_{c_{2}}+n_{v_{2}})n_{v_{2}}^{2}n_{c_{2}}^{0.5}\big{)}$.

Algorithm 3 AO Algorithm for Problem P2

1: Initialize a feasible initial point,
$(\boldsymbol{\eta}^{\mathrm{Sen}})^{(0)}$ and
$(\boldsymbol{\eta}^{\mathrm{Sen}})^{(0)}$.

2: Set the iteration number $n=1$.

3: repeat

4: Determine $(\boldsymbol{\eta}^{\mathrm{Com}})^{(n)}$ by using Algorithm 2.

5: Compute $(\boldsymbol{\eta}^{\mathrm{Sen}})^{(n)}$ by solving (23).

6: Set $n=n+1$.

7: until some stopping criterion is satisfied.

##  IV Numerical Results

We assume that the $M$ APs and $K_{d}$ users are uniformly distributed at
random within a square of size $D\times D~{}\text{km}^{2}$, whose edges are
wrapped around to avoid the boundary effects. The large-scale fading
coefficient $\beta_{mk}$ models the path loss and shadow fading, according to
$\beta_{mk}=\mathrm{PL}_{mk}10^{\frac{\sigma_{sh}z_{mk}}{10}}$, where
$\mathrm{PL}_{mk}$ represents the path loss, and
$10^{\frac{\sigma_{sh}z_{mk}}{10}}$ represents the shadow fading with the
standard deviation $\sigma_{sh}$, and $z_{mk}\sim\mathcal{CN}(0,1)$. We use
the three-slope model for the path-loss $\mathrm{PL}_{mk}$ (in dB) as

$\displaystyle\mathrm{PL}_{mk}\\!=\\!\\!\left\\{\begin{array}[]{ll}\\!\\!\\!\\!\\!-L\\!-35\log_{10}(d_{mk})&d_{mk}>d_{1},\\\
\\!\\!\\!\\!\\!-L\\!-15\log_{10}(d_{1})\\!-20\log_{10}(d_{mk})&d_{0}<d_{mk}\leq
d_{1}\\\
\\!\\!\\!\\!\\!-L\\!-15\log_{10}(d_{1})\\!-20\log_{10}(d_{0})&d_{mk}\leq
d_{0}.\end{array}\right.$

where $L$ is a constant depending on the carrier frequency, the user and AP
heights, given in [9]. We further use the correlated shadowing model for
$d_{mk}>d_{1}$ [9]. Here, we choose $\sigma_{sh}=8$ dB, $D=0.5$ km, $d_{1}=50$
m, and $d_{0}=10$ m. We further set the noise power $\sigma^{2}_{n}=-108$ dBm.
Let $\tilde{\rho}=1$ W and $\tilde{\rho}_{t}=0.25$ W be the maximum transmit
power of the APs and uplink training pilots, respectively. The normalized
maximum transmit powers ${\rho}$ and ${\rho}_{t}$ are calculated by dividing
these powers by the noise power $\sigma^{2}_{n}$.

In the absence of power control, termed as no power control (NPC) design, both
C-APs and S-APs transmit at full power. With NPC, the power coefficients at
the $m$-th C-AP are the same and
$\eta_{mk}=\big{(}N\sum_{k\in\mathcal{K}_{d}}\gamma_{mk}\big{)}^{-1}$,
$\forall k\in\mathcal{K}_{d}$. Moreover $\eta_{m}=1$, $\forall
m\in\mathcal{M}$. For comparison, two benchmark system designs are studied: 1)
Random AP mode selection with NPC (RAP-NPC), and 2) Greedy AP mode selection
with NPC (GAP-NPC).

Figure 2: The CDF of the minimum per user SE for different values of $\kappa$
($M=80$, $N=3$, $K_{d}=5$).

Figure 2 shows the cumulative distribution function (CDF) of the minimum per-
user SE for two different values of $\kappa$. For a fair comparison, the
achievable SE of the RAP-NPC is set to zero, when the MASR target value is not
satisfied. It is noteworthy that the $95\%$-likely minimum per-user SE
increases, when $\kappa$ decreases from $15$ to $10$. This behavior can be
explained by the fact that a great number of APs are assigned for downlink
communication by reducing $\kappa$. Furthermore, by applying power control,
the $95\%$-likely SE is improved by $50\%$ (for $\kappa=15$) and $35\%$ (for
$\kappa=10$) compared to the GAP-NPC design.

Figure 3: Average minimum SE versus MASR level ($K_{d}=5$, $NM=240$).

Figure 3 shows the minimum SE of the downlink communication system versus the
MASR requirement of the radar system, $\kappa$. In this figure, we assume that
a fixed number of service antennas, i.e., $MN=240$, is utilized to support
both communication and sensing applications. We observe that by increasing
$\kappa$, the achievable downlink SE is decreased since the number of C-APs is
decreased. Moreover, while RAP-NPC nearly fails to satisfy the MASR
requirements, the GAP scheme with NPC and OPC not only meets the MASR
requirements of the sensing operation, but also provides significant SE for
the communication system, especially when the number of antennas per each AP
is increased from $N=3$ to $N=12$.

##  V Conclusion

In this paper, we studied the problem of AP operation mode selection and power
control design in a cell-free massive MIMO ISAC system, aiming to support
multiuser downlink communication and single-target detection. After deriving
closed-form expressions for the per-user effective SINR and sensing MASR, we
formulated a max-min SE fairness problem. To deal with the complicated non-
convex problem, a greedy algorithm, grounded on long-term statistics, was
proposed for the AP operation mode design and an AO algorithm was developed
for power control design at all APs. Our numerical results highlighted that
our proposed GAP-OPC design significantly enhances the downlink SE compared to
the GAP-NPC and RAP-NPC benchmarks, while maintaining reliable sensing
performance.

##  Appendix A Proof of Proposition 2

In order to apply the use-then-forget technique to derive the downlink SE at
user $k$, we rewrite (II-D) as

$\displaystyle y_{k}$
$\displaystyle=\mathrm{DS}_{k}x_{c,k}\\!+\\!\mathrm{BU}_{k}x_{c,k}\\!+\\!\sum_{k^{\prime}\in\mathcal{K}_{d}\setminus
k}\\!\\!\\!\mathrm{IUI}_{kk^{\prime}}x_{c,k^{\prime}}\\!+\\!\mathrm{IR}_{k}x_{r}\\!+\\!n_{k},$
(24)

where

$\displaystyle\mathrm{DS}_{k}$
$\displaystyle\overset{\Delta}{=}\sum_{m\in\mathcal{M}}a_{m}\sqrt{\rho\eta_{mk}}\mathbb{E}\Big{\\{}\mathbf{g}^{T}_{mk}{\bf
t}_{mk}^{\mathrm{Com}}\Big{\\}}$ (25a) $\displaystyle\mathrm{BU}_{k}$
$\displaystyle\overset{\Delta}{=}\sum_{m\in\mathcal{M}}a_{m}\sqrt{\rho\eta_{mk}}\Big{(}{\bf
g}_{mk}^{T}{\bf t}_{k}-\mathbb{E}\big{\\{}{\bf g}_{mk}^{T}{\bf
t}_{k}\big{\\}}\Big{)},$ (25b) $\displaystyle\mathrm{IUI}_{kk^{\prime}}$
$\displaystyle\overset{\Delta}{=}\sum_{m\in\mathcal{M}}\sum_{k^{\prime}\in\mathcal{K}_{d}\setminus
k}a_{m}\sqrt{\rho\eta_{mk^{\prime}}}\mathbf{g}^{T}_{mk}{\bf
t}_{mk^{\prime}}^{\mathrm{Com}},$ (25c) $\displaystyle\mathrm{IR}_{k}$
$\displaystyle\overset{\Delta}{=}\sum_{m\in\mathcal{M}}(1-a_{m})\sqrt{\eta_{m}\rho}{\bf
g}_{mk}^{T}{\bf t}_{m}^{\mathrm{Sen}},~{}$ (25d)

respectively represent the strength of the desired signal ($\mathrm{DS}_{k}$),
the beamforming gain uncertainty ($\mathrm{BU}_{k}$), interference caused by
the $k^{\prime}$-th user ($\mathrm{IUI}_{kk^{\prime}}$) and the interference
caused by S-APs ($\mathrm{IR}_{k}$), respectively. By invoking (24), the
achievable downlink SE at the $k$-th user can be expressed as
$\mathrm{SE}_{k}=\Big{(}1-\frac{\tau_{\mathrm{p}}}{\tau}\Big{)}\log_{2}\big{(}1+\mathrm{SINR}_{k}\big{)}$,
where

$\displaystyle~{}\mathrm{SINR}_{k}\\!=\\!\frac{\big{|}\mathrm{DS}_{k}\big{|}^{2}}{\mathbb{E}\Big{\\{}\big{|}\mathrm{BU}_{k}\big{|}^{2}\Big{\\}}\\!+\\!\sum_{k^{\prime}\in\mathcal{K}_{d}\setminus
k}\mathbb{E}\Big{\\{}\big{|}\mathrm{IUI}_{kk^{\prime}}\big{|}^{2}\Big{\\}}\\!+\\!\mathbb{E}\Big{\\{}\big{|}\mathrm{IR}_{k}\big{|}^{2}\Big{\\}}\\!+\\!1}\\!.$
(26)

Now, we proceed to the desired signal term as

$\displaystyle~{}\mathrm{DS}_{k}$
$\displaystyle=\sqrt{\rho}\mathbb{E}\Big{\\{}\sum_{m\in\mathcal{M}}a_{m}\eta_{mk}^{1/2}(\hat{{\bf
g}}_{mk}+\tilde{{\bf g}}_{mk})^{T}\hat{{\bf g}}_{mk}^{*}\Big{\\}}$
$\displaystyle=\sqrt{\rho}\sum_{m\in\mathcal{M}}a_{m}N\eta_{mk}^{1/2}\gamma_{mk},$
(27)

where we have used the fact that $\hat{{\bf g}}_{mk}$ and $\tilde{{\bf
g}}_{mk}$ are zero mean and independent.

Noticing that the variance of a sum of independent RVs is equal to the sum of
the variances, we can derive
$\mathbb{E}\Big{\\{}\big{|}\mathrm{BU}_{k}\big{|}^{2}\Big{\\}}$ as

$\displaystyle\mathbb{E}\Big{\\{}\big{|}\mathrm{BU}_{k}\big{|}^{2}\Big{\\}}$
$\displaystyle=\rho\sum_{m\in\mathcal{M}}a_{m}\eta_{mk}\mathbb{E}\Big{\\{}\Big{|}{\bf
g}_{mk}^{T}\hat{{\bf g}}_{mk}^{*}-\mathbb{E}\Big{\\{}{\bf g}_{mk}^{T}\hat{{\bf
g}}_{mk}^{*}\Big{\\}}\Big{|}^{2}\Big{\\}}$
$\displaystyle=\rho\\!\\!\sum_{m\in\mathcal{M}}\\!\\!a_{m}\eta_{mk}\Big{(}\mathbb{E}\Big{\\{}\big{|}{\bf
g}_{mk}^{T}\hat{{\bf
g}}_{mk}^{*}\big{|}^{2}\Big{\\}}\\!-\\!\Big{|}\mathbb{E}\Big{\\{}{\bf
g}_{mk}^{T}\hat{{\bf g}}_{mk}^{*}\Big{\\}}\Big{|}^{2}\Big{\\}}\Big{)}$
$\displaystyle=\rho\\!\\!\sum_{m\in\mathcal{M}}\\!\\!a_{m}\eta_{mk}\Big{(}\mathbb{E}\Big{\\{}\big{|}\tilde{{\bf
g}}_{mk}^{T}\hat{{\bf
g}}_{mk}^{*}|^{2}\Big{\\}}\\!+\\!\mathbb{E}\Big{\\{}||\hat{{\bf
g}}_{mk}||^{4}\Big{\\}}\\!-\\!N^{2}\gamma_{mk}^{2}\Big{)}.$

By using the fact that $\mathbb{E}\Big{\\{}\big{|}\tilde{{\bf
g}}_{mk}^{T}\hat{{\bf
g}}_{mk}^{*}\big{|}^{2}\Big{\\}}=\mathbb{E}\Big{\\{}\tilde{{\bf
g}}_{mk}^{T}\mathbb{E}\big{\\{}\hat{{\bf g}}_{mk}^{*}\hat{{\bf
g}}_{mk}^{T}\big{\\}}\tilde{{\bf
g}}_{mk}^{*}\Big{\\}}=N\gamma_{mk}(\beta_{mk}-\gamma_{mk})$ and
$\mathbb{E}\Big{\\{}||\hat{{\bf
g}}_{mk}||^{4}\Big{\\}}=N(N+1)(\gamma_{mk})^{2}$, we get

$\displaystyle~{}\mathbb{E}\Big{\\{}\big{|}\mathrm{BU}_{k}\big{|}^{2}\Big{\\}}=\rho
N\\!\sum_{m\in\mathcal{M}}\\!a_{m}\eta_{mk}\gamma_{mk}\beta_{mk}.$ (28)

By following the same steps, we can obtain

$\displaystyle~{}\mathbb{E}\Big{\\{}\big{|}\mathrm{IUI}_{kk^{\prime}}\big{|}^{2}\Big{\\}}=\rho
N\sum_{m\in\mathcal{M}}a_{m}\eta_{mk^{\prime}}\gamma_{mk^{\prime}}\beta_{mk}.$
(29)

Moreover, by substituting (7) into (25d) and noticing that ${\bf
t}_{m}^{\mathrm{Sen}}({\bf t}_{m}^{\mathrm{Sen}})^{H}=\frac{1}{N}{\bf I}_{N}$,
we get

$\displaystyle~{}\mathbb{E}\Big{\\{}\big{|}\mathrm{IR}_{k}\big{|}^{2}\Big{\\}}$
$\displaystyle=\sum_{m\in\mathcal{M}}(1-a_{m})\mathbb{E}\Big{\\{}|{\bf
g}_{mk}^{T}{\bf t}_{m}^{\mathrm{Sen}}|^{2}\Big{\\}}$
$\displaystyle=\rho\sum_{m\in\mathcal{M}}{\eta_{m}}(1-a_{m})\beta_{mk}.$ (30)

To this end, by substituting (A), (28), (29), and (A) into (26), the desired
result in (13) is obtained.

## References

  * [1] F. Liu, Y. Cui, C. Masouros, J. Xu, T. X. Han, Y. C. Eldar, and S. Buzzi, “Integrated sensing and communications: Toward dual-functional wireless networks for 6G and beyond,” IEEE J. Sel. Areas Commun., vol. 40, no. 6, pp. 1728-1767, Jun. 2022.
  * [2] A. Liu et al., “A survey on fundamental limits of integrated sensing and communication,” IEEE Commun. Surveys Tuts., vol. 24, no. 2, pp. 994- 1034, Second Quarter 2022
  * [3] F. Liu, C. Masouros, A. Li, and T. Ratnarajah, “Robust MIMO beamforming for cellular and radar coexistence,” IEEE Wireless Commun. Lett., vol. 6, no. 3, pp. 374-377, Jun. 2017.
  * [4] B. Li and A. P. Petropulu, “Joint transmit designs for coexistence of mimo wireless communications and sparse sensing radars in clutter,” IEEE Trans. Aerosp. Electron. Syst., vol. 53, no. 6, pp. 2846-2864, Dec. 2017
  * [5] M. Elfiatoure, H. Q. Ngo, and M. Matthaiou, “Coexistence between massive MIMO and radar communications: Performance analysis,” J. Commun. Inf. Networks, vol. 8, no. 1, pp. 37-47, Mar. 2023.
  * [6] A. R. Chiriyath, B. Paul, and D. W. Bliss, “Radar-communications convergence: Coexistence, cooperation, and co-design,” IEEE Trans. Cognitive Commun. and Networking, vol. 3, no. 1, pp. 1-12, Mar. 2017.
  * [7] F. Liu, C. Masouros, A. Li, H. Sun, and L. Hanzo, “MU-MIMO communications with MIMO radar: From co-existence to joint transmission,” IEEE Trans. Wireless Commun., vol. 17, no. 4, pp. 2755-2770, Apr. 2018.
  * [8] M. Matthaiou et al., “The road to 6G: Ten physical layer challenges for communications engineers,” IEEE Commun. Mag., vol. 59, no. 1, pp. 64– 69, Jan. 2021.
  * [9] H. Q. Ngo, A. Ashikhmin, H. Yang, E. G. Larsson, and T. L. Marzetta, “Cell-free massive MIMO versus small cells,” _IEEE Trans. Wireless Commun.,_ vol. 16, no. 3, pp. 1834-1850, Mar. 2017.
  * [10] Z. Behdad, O. T. Demir, K. W. Sung, E. Björnson, and C. Cavdar, ¨ “Power allocation for joint communication and sensing in cell-free massive MIMO,” in Proc. IEEE GLOBECOM, Dec. 2022, pp. 4081-4086.
  * [11] U. Demirhan and A. Alkhateeb, “Cell-free ISAC MIMO systems: Joint sensing and communication beamforming,” arXiv preprint arXiv:2301.11328, 2023\.
  * [12] T. L. Marzetta, “Noncooperative cellular wireless with unlimited numbers of base station antennas,” IEEE Trans. Wireless Commun., vol. 9, no. 11, pp. 3590-3600, Nov. 2010.
  * [13] J. A. Sutton, H. Q. Ngo, and M. Matthaiou, “Hardening the channels by precoder design in massive MIMO with multiple-antenna users,” IEEE Trans. Veh. Technol., vol. 70, no. 5, pp. 4541-4556, May 2021.
  * [14] S. Boyd and L. Vandenberghe, Convex Optimization. New York, NY: Cambridge University Press, 2004.
  * [15] H. H. M. Tam, H. D. Tuan, D. T. Ngo, T. Q. Duong, and H. V. Poor, “Joint load balancing and interference management for small-cell heterogeneous networks with limited backhaul capacity,” IEEE Trans. Wireless Commun., vol. 16, no. 2, pp. 872-884, Feb. 2017.