# GaitSADA: Self-Aligned Domain Adaptation for mmWave Gait Recognition

Ekkasit Pinyoanuntapong1, Ayman Ali1, Kalvik Jakkala1, Pu Wang1, Minwoo Lee1,
Qucheng Peng2, Chen Chen2, Zhi Sun3  1University of North Carolina at
Charlotte
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>2University of Central Florida
<EMAIL_ADDRESS><EMAIL_ADDRESS>3Tsinghua University
<EMAIL_ADDRESS>
###### Abstract

mmWave radar-based gait recognition is a novel user identification method that
captures human gait biometrics from mmWave radar return signals. This
technology offers privacy protection and is resilient to weather and lighting
conditions. However, its generalization performance is yet unknown and limits
its practical deployment. To address this problem, in this paper, a non-
synthetic dataset is collected and analyzed to reveal the presence of spatial
and temporal domain shifts in mmWave gait biometric data, which significantly
impacts identification accuracy. To mitigate this issue, a novel self-aligned
domain adaptation method called GaitSADA is proposed. GaitSADA improves system
generalization performance by using a two-stage semi-supervised model training
approach. The first stage employs semi-supervised contrastive learning to
learn a compact gait representation from both source and target domain data,
aligning source-target domain distributions implicitly. The second stage uses
semi-supervised consistency training with centroid alignment to further close
source-target domain gap by pseudo-labelling the target-domain samples,
clustering together the samples belonging to the same class but from different
domains, and pushing the class centroid close to the weight vector of each
class. Experiments show that GaitSADA outperforms representative domain
adaptation methods with an improvement ranging from 15.41% to 26.32% on
average accuracy in low data regimes. Code and dataset will be available at
https://exitudio.github.io/GaitSADA

###### Index Terms:

mmWave radar, Domain Adaptation, Deep Learning, Gait Recognition

##  I Introduction

User identification is an essential component for many Internet of Things
(IoT) applications, e.g. personalized environmental control, security
management, and access control for automatic doors. Short-distance biometrics
such as fingerprints, face, iris, palm, and finger vein patterns typically
demand active user participation within a certain proximity. In contrast,
gait, i.e., human walking style, can be exploited for user identification from
a distance without the subject’s participation or interference. Moreover, it
is shown in existing study [1] that human gait is very hard to spoof or mimic
because our own gait works against us when we try to imitate someone elses
gait.

The majority of gait recognition techniques are vision-based, which directly
perform identification tasks on video data. Appearance-based and model-based
approaches are two general approaches for vision-based gait recognition.
Appearance-based approaches [2][3][4] utilize silhouettes from a video
sequence obtained from background subtraction. Although background subtraction
is accurate in a lab context, it is difficult to use in a cluttered and
dynamic real-world environment. Model-based approaches [5][6] try to tackle
this issue by employing the existing pose estimators to improve robustness to
environmental variation. However, all vision-based approaches require access
to RGB visual images which are prone to performance loss from variations of
ambient light settings and the social restriction from privacy concerns.

Radar-based gait recognition is an emerging user identification solution that
exploits radar return signals to capture human gait biometric [7, 8, 9, 10,
11, 12, 13]. Unlike camera surveillance systems, radar-based gait recognition
systems have several unique advantages. First, they are privacy-preserving as
they do not rely on images of human subjects. Next, they can operate in
adverse weather and lighting conditions such as heavy fog, smoke, rain, and
zero-light conditions. Finally, they can see through some opaque objects
depending on its material and the frequency of the radar [14, 15]. Radar-based
methods operate by emitting modulated electromagnetic signals, which are
scattered and returned back by the subjects in the signal propagation path.
The radar return signals capture the micro-Doppler (mD) signature, which can
be represented as radar spectrogram. mD signature or radar spectrogram
documents the time-dependent Doppler frequency shifts, which are caused by the
small-scale vibrations or rotations of human body parts during walking. Each
person has his/her unique walking style, i.e., gait, which leads to the unique
mD signature that can be exploited for gait recognition.

Existing radar-based gait recognition systems generally operate either in
sub-6GHz WiFi bands [7, 8, 9] or the emerging mmWave bands (30 to 300 GHz)
[10, 11, 12, 13]. Compared with sub-6GHz radar systems, mmWave radar can
promise much fine-grained motion capture with significantly enhanced velocity
and location resolutions. This feature is enabled by orderly shorter
wavelength and wider bandwidth of mmWave radar signals, compared with the
sub-6GHz counterparts. Besides promising motion capture capacities, mmWave
radar system only needs a single radar device that integrates both radar
transmitter and receiver, while sub-6GHz radar system is generally implemented
by a pair of WiFi devices to act as distributed radar transmitter and
receiver, which have to run modified WiFi firmware to extract mD signatures
from CSI signals.

The superior gait recognition performance of mmWave radar systems has been
demonstrated in our previous research [10] and other related work [11, 12,
13]. However, the spatial and temporal generalization performance of this
technology is still unknown. We found that deep neural networks (DNNs) trained
on radar-based gait data suffer from substantial performance degradation with
the progression of time and in new locations. Existing radar-based gait
recognition systems proposed in the literature, test their performance on a
dataset collected alongside the training data and fail to uncover the temporal
degradation. Furthermore, they seldom address the spatial generalization
issues from introducing new environments. The root cause of this
generalization issue is domain shift [16]. Simply put, it occurs when the
testing data distribution differs from that of the training data, i.e., the
i.i.d. assumption is violated. Domain shift problem can be further exacerbated
by the harmful environment reflections induced by mmWave signals. These
environment reflections include static ones, which directly bounce from
stationary obstructions e.g. walls, desks, and chairs, and dynamic ones, which
indirectly bounce from other stationary obstructions and then bounce from
human bodies. The harmful environment reflections carry delayed and distorted
mD signature information, which cannot be completely removed. This leads to
environment-induced temporal domain shift (TDS) and spatial domain shift
(SDS). Moreover, human gait, although unique to each individual, could contain
many variations, potentially as a consequence of a subject’s mood and clothing
changes, further aggravate TDS and SDS.

In this paper, we tackle the challenge of spatial-temporal domain drift
problem in mmWave radar-based gait recognition. To address this issue, we
propose a novel domain adaptation method, GaitSADA. GaitSADA is based on
simple distribution alignment components and consists of a two-stage approach:
pretraining and fine-tuning. In the pretraining stage, the model learns
salient gait representations and noise invariance from both source and target
domains using semi-supervised contrastive learning. In the fine-tuning stage,
the model generalizes the learned gait representation to the target domain by
leveraging pseudo-labels, semi-supervised consistency training, and centroid
alignment. Experiments demonstrate that GaitSADA outperforms classic
adversarial domain adaptation methods and achieves significant accuracy
improvements in low data regimes and challenging domain drift conditions. Our
contributions can be summarized as follows:

  * •

We curate a non-synthetic dataset consisting of mmWave radar-based gait
biometric data. This dataset, for the first time, allows one to study and
improve the spatio-temporal generalization performance of a radar-based
biometric identification system.

  * •

We perform extensive research on the SDS and TDS from four different location
environments on different days, experimenting with five state-of-the-art
(SOTA) domain adaptation methods to mitigate these drifts and show the issue
of accuracy drop in low data regimes of SOTA domain adaptation methods.

  * •

We propose a novel domain adaptation method for mmWave gait recognition
system, self-aligned domain adaptation for mmWave gait recognition (GaitSADA),
by jointly exploiting semi-supervised contrastive learning, semi-supervised
consistency training, centroid alignment and deliberately designed radar
spectrogram augmentations. GaitSADA outperforms existing SOTA methods in
1-day, 2-day, and 3-day data by 15.41%, 5.83%, and 2.63% respectively on
average of domain drift cases.

##  II Related Work

Unsupervised domain adaptation (UDA), semi-supervised learning, and self-
supervised learning are common strategies to reduce the dependency on
expensive manually annotated data. UDA aims to mitigate the domain gap by
transferring knowledge learned from label data in source domain and transfer
to unlabeled data in target domain, while semi-supervised learning employs a
few manually annotated samples and generalizes to large unlabeled samples
without domain shift assumption. Self-supervised utilizes the unlabeled data
to pre-train the model to learn the general semantics of the data.

Unsupervised domain adaptation can tackle the domain drift and the lack of
labeled data problem by minimizing the discrepancy between source and target
distributions without using expensive labeled data from the target domain.
Many UDA techniques, inspired by Generative Adversarial Networks (GAN), apply
an adversarial loss to confuse the domain discriminator network that tries to
recognize the domain difference. ADDA [17] exploits two networks for source
and target domains, respectively and forces the output representation to be in
the same distribution. Domain-adversarial neural network (DANN) trains domain
discriminator along with the feature extractor that confuses the domain
discriminator by integrating a gradient reversal layer (GRL) [18]. Although
the overall distribution is aligned, the class-wise discrepancy can still
cause performance drop in the classification tasks. Conditional Adversarial
Domain Adaptation (CDAN) [19] takes this into account by implementing
multilinear conditioning of the feature representation and the class-wise
output prediction.

Semi-supervised methods decrease the dependency on labeled data by combining
normally smaller sets of labeled data with the vast amounts of unlabeled data
that are available in many use cases. Since the classification loss is unknown
for unlabeled data, the state-of-the-art (SOTA) semi-supervised Learning
solutions, such as UDA [20], FixMatch [21], MixMatch [22], and RemixMatch [23]
mainly rely on consistency training, which learns target unlabeled data
$u_{i}$ by generating pseudo-labels from models’ predictions on weakly-
augmented data and filtering only high confidence data by selecting the
prediction that has confidence more than the threshold. The cross-entropy $H$
is applied against the pseudo labels and the model’s prediction of the highly-
augmented versions of the same unlabeled data.

Similar to semi-supervised learning, self-supervised Learning utilizes large
unlabeled data to improve performance from label data trained in a supervised
learning fashion. However, self-supervised learning normally pre-trains
unlabeled data to learn general data semantics before transferring pre-trained
weights to further learn label data in a supervised learning fashion in the
next stage. Many research [24], [25], [26], [27] enforce contrastive loss to
learn the structure of the data without labels by imposing different
augmentation to the same data and learning the similarity between the features
of same data (positive samples). The challenge is that such similarity loss
functions can trigger a dimensional collapse in which the model outputs the
exact same embedding for all input. While MoCo [27] and BYOL [25] incorporate
momentum encoder and stop gradient to prevent the collapse, SimCRL [24] and
CPC [26] employ negative samples to discriminate the features of different
data. SimCLR [24] proposed NT-Xent (the normalized temperature-scaled cross-
entropy loss) which simply maximizes the similarity between the feature of the
same data (positive samples) and minimizes the different ones (negative
samples).

##  III mmWave gait biometric Dataset

The dataset is created by exploiting the mmWave radar-based gait recognition
system, STAR, developed in our previous work [10]. In this section, the
preliminaries of mmWave radar sensing are first presented, the data collection
and gait recognition procedure using STAR system is then introduced, and
finally spatial-temporal domain drift issues are revealed.

Figure 1: Overview of STAR system

###  III-A STAR System Overview

STAR is the first one in the literature to exploit deep learning for mmWave
radar-based gait recognition. As shown in Fig. 1, STAR consists of three
components: radar spectrogram acquisition, spectrogram enhancement and deep
learning(DL)-based gait spectrogram recognition. The spectrogram acquisition
module tracks and captures the radar signal returns that directly bounce off
from human body, while filtering out undesired environment-induced radar
returns. This is achieved by applying a sequence of radar widely-adopted radar
signal processing algorithms including one-dimensional Fourier transform (ID-
FFT) for range estimation, the mean subtraction for static reflection removal,
CA-CFAR for robust target detection against a background of noise, clutter and
interference, Doppler estimation via 2D-FFT to estimate target velocity, and
sliding-window based tracking method to predict the target future state (e.g.,
location and velocity). The target tracking module acts a moving spatial
passband filter so that only the radar returns from the target can be
received. Performing short-time Fourier transform (STFT) on received radar
returns generates the raw radar spectrogram, whose quality is further improved
by exploiting the spectrogram enhancement module that consists of a sequence
of frequency-domain noise removal processes [10]. The enhanced spectrograms
are then encoded by a CNN-based feature encoder that maps the high-dimension
spectrogram data into the compact low-dimension feature vector. Finally, the
compact feature vectors are used to train a multiclass classifier that
predicts the corresponding user identities.

###  III-B Data Collection Setup

We collected gait data from 10 volunteers between the ages of 18-35. Each
subject’s data was collected in four different locations. The source location
was a research space with cubicles, and the other three areas consisted of a
server, conference, and an office room. By maintaining four distinct
locations, we introduce SDS. In the source location, data was collected on 10
different days for each subject. 5 separate days of data was acquired for each
of the three other locations, which are used as target domains. A participant
can either walk towards the radar or walk away from the radar, each of which
is counted as one walking instance and generates one spectrogram data sample.
The data collection was limited to 100 data samples per person in the source
location and 50 data samples per person in each of the target locations on any
given day. The gait spectrogram samples of a person are shown in Fig. 2.

###  III-C Establishing the Presence of TDS and SDS

We study the presence of TDS and SDS issues by training gait recognition model
using the labelled data from source domain and test the model performance
across different target domains. In particular, the feature encoder is a
modified version of ResNet-50 [28], as described in section V-B, to extract
features from mmWave biometric data samples. The loss function is a cross-
entropy loss for multiclass classification. To establish the presence of SDS,
we train the model on the data from 1 to 3 days of source location and test
the model on the data from other three unseen locations. SDS is evident from
the performance degradation in all target domains as shown in table I. To
reveal the existence of TDS, we train the model on the data from 1 to 3 days
of the source location and test the model on the data from the remaining days
of the source location. TDS can be observed due to the degraded testing
accuracy.

Figure 2: Gait spectrograms of the same person for 5 days at different
locations, source (laboratory), conference room, server room, and office.
TABLE I: The result of supervised learning trained on source data and tested
on different days of data (i.e., temporal) and other different environments
(i.e., server, conference, office). The accuracy drop reveals the domain
discrepancy when the experiments are conducted on different days and
environments. The results report in % of accuracy.

# of training days (Source) Temporal Server Conference Office 1-day 76.05
34.93 47.00 33.39 2-days 88.46 53.05 71.00 49.92 3-days 95.25 71.77 85.20
67.61

##  IV Method: GaitSADA

###  IV-A Method Overview

Our novel method, self-aligned domain adaptation for mmWave gait recognition
(GaitSADA), jointly exploits self- and semi-supervised learning techniques to
realize the unsupervised domain adaptation by promoting intra-class
compactness and inter-class distance while generalizing knowledge from label
data to unlabeled data using pseudo-label. In particular, we leverage the
unlaballed data from the target domains in two training stages: pretaining and
fine-tuning as illustrated in Fig. 3.

During pretraining stage, we propose a semi-supervised contrastive learning
paradigm to learn the salient gait representation from a mixed dataset that
consists of data from both source and and target domains. Our semi-supervised
contrastive learning can explicitly avoid collapse problem without requiring
the sophisticated techniques adopted by classic self-supervised (unsupervised)
contrastive learning paradigms, e.g., negative samples, stop gradients, and
moment encoders. This is achieved by maximizing the agreement between
representations produced by encoders fed with different augmentations of the
same spectrogam in the mixed dataset, while enforcing the encoder fed with
labelled data in source domain to produce discriminative representations with
correct class predictions.

After the pretrained stage, the feature encoder implicitly aligns the data
distributions of source domain and target domain in the feature space. As a
result, the representations produced by such encoder are discriminative enough
to provide relatively high-confidence pseudo labels for unlabelled data in the
target domains, which is essential for the semi-supervised consistency
training in the fine-tuning stage that can effectively minimize the embedding
vectors of the target domain data from the same class to be close to each
other. However, classic semi-supervised consistency training only shows the
best performance when the label and unlabelled data come from the same
distribution. To address this issue, we propose the centroid alignment scheme
that can be integrated with semi-supervised consistency training to mitigate
source-target domain gap.

###  IV-B Pretraining via Semi-supervised Contrastive Learning

Supervised Learning with Additive Margin Softmax The main goal of this
component is to promote intra-class compactness and inter-class separability.
As the regular softmax loss is effective only in maximizing the inter-class
difference, it is ineffective at minimizing the intra-class variation.
Therefore, we adopt Additive Margin Softmax loss function [29] to bring the
classification boundary closer to the weight vector of each class by adding a
margin $m$ to the boundary as

$L_{supervised}=-\frac{1}{B}\sum_{i=1}^{B}log\frac{e^{s\cdot(cos\theta_{y_{i}})-m}}{e^{s\cdot(cos\theta_{y_{i}})-m}+\sum_{j=1,j\neq
y_{i}}^{C}e^{s\cdot cos\theta_{j}}},$ (1)

where
$cos\theta_{y_{i}}=\frac{W^{T}_{y_{i}}f_{i}}{\left\|W_{y_{i}}\right\|\left\|f_{i}\right\|}$
represents cosine similarity of the inner product between the normalized
weight $W_{y_{i}}$ and the normalized feature $f_{i}$ of class $i$. Additional
parameter $s$ is a scale parameter to control the temperature.

Figure 3: Overall Architecture of GaitSADA. (a) Stage 1. Pretraining via semi-
supervised contrastive learning to learn a compact gait representation from
both source and domain data, which also implicitly mitigates the source-domain
distribution shift. (b) Stage 2. fine-tuning via semi-supervised consistency
training with centroid alignment. The fine-tuning stage aims to further
improve model generalization performance by pseudo-labelling the target-domain
samples, clustering together the samples belonging to the same class but from
different domains, and pushing the class centroid close to the weight vector
of each class. The encoders are sharing the weights.

Positive Unlabelled Contrastive Learning. The objective of the second
component is to learn noise invariance by training the model to generate the
same features regardless of the artificial noise that is added to the data
without relying on label data. We carefully design multiple augmentations to
simulate the radio wave noise such as frequency noise, temporal noise,
different resolutions, and white noise as described in Section V-A. The model
learns noise invariance by minimizing cosine similarity between features from
raw data $f$ and its augmented version $\hat{f}$ with the loss

$L_{similarity}=\frac{1}{B}\sum_{i=1}^{B}\frac{f_{i}^{T}\hat{f_{i}}}{\left\|f_{i}\right\|\left\|\hat{f_{i}}\right\|}.$
(2)

This loss is inspired by NT-Xent lose from SimCLR [24] which is the
contrastive self-supervised learning method that minimizes cosine similarity
of positive samples (i.e., different augmentation of the same data) and
maximizes the cosine similarity of negative samples (i.e., different data).
The objective of negative samples in NT-Xent is to prevent dimensional
collapse. However, since our self alignment loss is trained jointly with
supervised learning loss, this can explicitly avoid the collapse issue without
relying on negative samples. This is because if the encoder produces constant
embedding vectors for different spectrograms (i.e., dimensional collapse
occurs) via self-supervised learning branch, this will prevent supervised
learning from correctly predicting the classes of the labelled spectrograms.
As a result, we are able to retain only the positive similarity and eliminate
the negative marginal fraction. The final loss function in the first stage
combines two loss functions to encourage the intra-class concentration and
inter-class distance while being invariant to perturbation from noise by
utilizing unlabeled data from target domain:

$L_{stage1}=L_{supervised}+L_{similarity}.$ (3)

###  IV-C Fine-tuning via Consistency Training with Centroid Alignment

In the first stage, the model has learned the well-aligned class-wise feature
from label data in the source domain and unlabelled data from target domain.
The objective of the second stage is to further regulate the feature of the
target domain data of the same class to be close to each other.

Consistency Training. Unlike the source domain, the labels in the target
domain are unknown. Therefore, we cannot directly align target feature by
class. To address this problem, we enforce consistency training via pseudo
labels instead of the true labels of the target domain. Consistency training
utilizes unlabeled data by relying on the assumption that when given perturbed
variations of the same input, the model ought to provide similar predictions.
The two identical models are used for consistency training. The first branch
takes raw input data with no augmentation in order to provide accurate pseudo
labels. We compute the model’s predicted class distribution of a given
unlabeled spectrogram input: $q_{i}=p(y_{i}|u_{i})$ as soft pseudo labels. The
cross-entropy $H$ is applied against the pseudo labels distribution and the
model’s prediction distribution of the augmented unlabeled data
$\mathcal{A}(u_{i})$. The loss is defined as

$L_{consistency}=\frac{1}{B}\sum_{i=1}^{B}\mathds{1}(\arg\max_{y}(q_{i})\geq\tau)H(q_{i},\hat{q}_{i}),$
(4)

where $\hat{q}_{i}=p(y_{i}|\mathcal{A}(u_{i}))$ represents the probability
distribution of class $y_{i}$ given augmented input $\mathcal{A}(u_{i})$ and
$\tau$ is a scalar parameter denoting the threshold to filter only the high
confidence prediction. To maintain the same standard of Additive Margin
Softmax (Eq. (1)) in the first stage, we normalize weights and features of the
classifier which can be referred to as cosine similarity $cos\theta$. The
scale $s$ remains the same, but we remove margin $m$ for accurate pseudo
labels:

$p(y_{i}|\cdot)=\frac{e^{s\cdot(cos\theta_{y_{i}})}}{\sum_{j=1}^{c}e^{s\cdot
cos\theta_{j}}},$ (5)

where $c$ is the number of classes. The predictions are likely to be
inaccurate at the beginning of the training. For this, most semi-supervised
learning algorithms [30][31][22][23] initialize training with a small weight
$\lambda$ for the consistency loss term and increase it over time. However,
with the confidence threshold, the pseudo label produces a natural curriculum
learning [21]. As less pseudo-label data are available at first since the low
confidence predictions are filtered out. When the model is more accurate, it
then provides more predictions with a high degree of confidence. As a result,
the threshold automatically increases the number of pseudo-label data without
explicitly specifying incremental weight $\lambda$.

Centroid Alignment. While consistency training loss tends to move weights $W$
of the classifier toward target features, supervised learning loss still
attempts to bring the weights back to source features. The weights can
oscillate between these two groups as a result of the domain gap between the
source and target features as shown in Fig. 4. To mitigate this issue, we
impose centroid alignment loss to gradually bring the weights to the center of
source and target features by applying Additive Margin Softmax to class-wise
weights and class centroids as

Figure 4: Centroid Alignment. The centroid of source features tends to be
closed to weight vector $W$ of the corresponding class following the
supervised learning loss. On the other hand, the centroid of target features
is likely to be further away due to the domain gap. Centroid Alignment
encourages weight vector $W$ to move to a new centroid which is the center of
both source and target features that can be a better representation as more
target features are learned from consistency training.

$L_{centroid}=-\frac{1}{C}\sum_{c=1}^{C}log\frac{e^{s\cdot(cos\theta_{c})-m}}{e^{s\cdot(cos\theta_{c})-m}+\sum_{j=1,j\neq
i}^{C}e^{s\cdot cos\theta_{j}}},$ (6)

where
$cos\theta_{c}=\frac{W^{T}_{c}Z_{c}}{\left\|W_{c}\right\|\left\|Z_{c}\right\|}$
is the inner product of normalized weights $W_{i}$ and centroids $Z_{i}$ of
each class. The class-wise centroid is from the Exponential Moving Average
(EMA) of both source and target features with hyperparameter $\alpha$ to
control the rate of change as

$Z_{c_{i}}=\alpha\hat{Z}_{c_{i}}+(1-\alpha)Z_{c_{i-1}}.$ (7)

The current batch centroid $\hat{Z}_{c_{i}}$ is the average features of source
$f_{s_{i}}$ and target $f_{t_{i}}$. However, because the true class of the
target domain is unknown, we cannot directly assign features to the class
centroid. Instead, the hardmax (i.e., argmax) of the model’s prediction is
used as a pseudo label of the class of each sample. Also, target features can
be noisy at the beginning of the training, similar to consistency loss (Eq.
(4)), only the high confidence labels are used for calculating the centroids
to avoid false predictions of pseudo labels. Therefore,

$\hat{Z}_{c_{i}}=\frac{\sum_{i=1}^{B}f_{s_{i}}+\mathds{1}(\arg\max_{y}(q_{i})\geq\tau)f_{t_{i}}}{B+\sum_{i=1}^{B}\mathds{1}(\arg\max_{y}(q_{i})\geq\tau)}.$
(8)

Lastly, we impose Additive Margin Softmax for source label data as a
supervised learning loss along with consistency training and Centroid
Alignment loss to prevent dimensional collapse. So the final loss equation is
shown as

$L_{stage2}=L_{supervised}+L_{consistency}+\lambda L_{centroid}.$ (9)

##  V Experiment Setup

Texas Instruments (TI) IWR1642EVM boost board interfaced with a DCA1000EVM
board to collect the raw mmWave radar signal. The radar system consists of two
transmitting antennas, four receiving antennas, and 120° view of the azimuth
plane. The radar system supports up to a 4Ghz bandwidth operating on 77 GHz to
81 GHz. To configure FMCW wave parameters such as chirp width, repetition
time, and chirp slope from our radar device, we use a Dell Latitude 7480
laptop with TI mmWave studio software as a control system.

Laboratory location data is used as source data and the other location data
(i.e. conference room, server room, and office) is referred to as target data.
For the spatial domain drift experiment, day 1 to day 3 samples are used as a
training set for all source and target locations. The rest of the data is used
for testing. For the temporal domain drift experiment, we employ 1 to 3 days
of laboratory location as source data and the next consecutive days of
laboratory location as target data. For example, considering the temporal
2-day case, the first and the second-day data of laboratory location is a
source domain, the third and the fourth-day data is a target domain, and the
fifth to the tenth-day data is utilized as a test set.

###  V-A Spectrogram Augmentation

(a) Cutout
Horizontal
.

(b) Cutout
Vertical
.

(c) White Noise
.
.

(d) Zoom
.

(e) Shear
.

(f) Rotation
.

Figure 5: Spectrogram augmentations. TABLE II: Result of training on 1 to 3
days on the data from laboratory location (source domain) while adapting to
different 1 to 3 days data of same location (i.e., temporal target domain) and
1 to 3 days of different target locations, (i.e., server, conference, and
office) with 5 representative methods, GRL [18], GAN [32], ADDA [17], CDAN
[19], and FixMatch [21] along with ”Supervised Learning” which train only on
source label data without domain adaptation and lastly, our proposed method,
GaitSADA. The results report in % of accuracy.

# of training days Methods Temporal Server Conference Office Average 1-day
Supervised 76.05 34.93 47.00 33.39 47.84 GAN 79.57 38.14 66.00 32.93 54.16 GRL
74.97 40.14 62.10 36.14 53.34 ADDA 80.96 37.71 65.21 32.60 54.12 CDAN 73.67
45.15 62.40 43.24 56.12 FixMatch 94.32 44.24 69.40 49.05 64.25 GaitSADA (our)
96.3 75.98 82.6 63.76 79.66 2-day Supervised 88.46 53.05 71.00 49.92 65.61 GAN
92.18 59.96 77.80 59.76 72.43 GRL 85.98 69.07 80.90 59.96 73.98 ADDA 88.55
51.04 75.63 55.62 67.71 CDAN 91.72 63.36 85.90 63.16 76.04 FixMatch 95.82
75.38 90.10 82.58 85.97 GaitSADA (our) 97.63 79.88 99.1 90.59 91.8 3-day
Supervised 95.25 71.77 85.20 67.61 79.96 GAN 96.79 78.18 94.30 85.89 88.79 GRL
96.11 84.58 89.30 75.58 86.39 ADDA 97.38 79.79 94.79 88.13 90.02 CDAN 98.03
91.39 93.80 86.19 92.35 FixMatch 99.29 93.29 98.30 88.89 94.94 GaitSADA (our)
99.68 93.69 99.2 97.7 97.57

Vertical noise stripe (Fig. 5(b)) simulates missing information in the random
time interval of the spectrogram. Similarly, the horizontal noise stripe (Fig.
5(b)) represents the missing information in some frequency range. We implement
horizontal and vertical noise together by randomly adding two of these noises
per data with a probability of 1/3 to be horizontal and 1/3 to be vertical,
the rest 1/3 is for no cutout noise. The small side of the stripe is random
between 2 to 8 pixels. The third augmentation is pixel-wised random numbers
adjusted to the whole data uniformly simulating white noise (Fig. 5(c)) in a
spectrogram. The uniformly random noise from -0.5 to 0.5 is added to the data
with a probability of occurrence 2/3. Next is the zoom-in and out effect.
Zoom-in will drop some information from the start and the end, also the lowest
and the highest frequency of the spectrogram. Zoom effects (Fig. 5(d))
represent the data in different resolutions. The scale is applied to the data
between 0.8 to 1.2 randomly. These four augmentations act as virtual noise to
simulate the noisy real-world data thereby improving domain generalization.
Not only these meaningful noises, but we also experiment with regular computer
vision augmentations and we found that shearing (Fig. 5(e)) and rotation (Fig.
5(f)) help boost performance as well. The shear angle in a counter-clockwise
direction in degrees is random from 0 to 5. Meanwhile, 5-degree ranges are
randomly applied for rotation. We apply these augmentations to our preprocess
pipeline and then apply all augmentation to every spectrogram data with random
magnitude. In addition to consistency training, augmentation is used in
supervised learning during the first training stage.

###  V-B Model and Training Details

A feature extraction model is modified from Resnet-50 to support one channel
instead of three channels and use SELU as an activation function. The other
parts are similar to the original Resnet-50 [28] composed of five stages of
convolution blocks and identical blocks including a stack of convolutional
blocks, batch normalization, and activation with skip connections in each
block and output 128 features. The classification model is simply a linear
layer with 10 neurons to classify 128 features from feature extraction into 10
classes.

We apply the same settings for all 1- to 3-day cases and all locations in both
stages 1 and 2 to be consistent in all domain configurations. Every
convolution layer applies $1\mathrm{e}{-4}$ of $l^{2}$ kernel regularizer and
trains for 10,000 epochs. The training optimizer is ADAM with an initial
learning rate starting at $1\mathrm{e}{-3}$ and ending at $1\mathrm{e}{-5}$,
applying the polynomial learning rate decay scheduler with learning rate
restart. We use a batch size of 64. The confidence threshold of the
consistency pseudo-label is 0.97. Additive Margin Softmax Loss applies
$l^{2}$-normalization to both input and weights to the fully connected layer
head without bias, with $s=10$ and $m=0.2$. All losses in both stages 1 and 2
are applied with equal weights except centroid alignment loss with $\alpha$ =
0.05 to gradually move centroids of the classes without breaking the current.

###  V-C Evaluation Results

We evaluate the model without any domain adaptation (supervised learning) as a
baseline comparing with the five representative domain adaptation methods
including GAN [33], GRL [18], ADDA [17], CDAN [19], and FixMatch [21] along
with our method, GaitSADA, on the four different locations (laboratory,
conference room, server room, and office) to explore the spatial domain
adaptation using 1-3 days of data for training and the rest of data for
testing. Meanwhile, We also study naive supervised learning without domain
adaptation to understand the domain discrepancy issue in low data regime,
temporal case, and spatial case.

As summarized in Table II, our GaitSADA method outperforms the other five
representative domain adaptation methods on all spatial domain adaptation
(i.e., different locations) and temporal domain adaptation (i.e., different
days on laboratory location). Without applying domain adaptation (i.e.,
supervised training in Table II), the performance of 3-day cases is 79.96% on
average for all locations and temporal case, then drop to 65.61% and 47.85% in
2-day and 1-day data which apparently shows that the amount of data helps
mitigate the domain shift issues in both spatial and temporal cases.

GAN, GRL, and ADDA yield similar results and help alleviate the domain drift
issue by boosting accuracy from the supervised method on average in all 3-day
cases. However, we observe the accuracies drop compared to the supervised
learning method in 1-day case for all these three methods, e.g., GAN and ADDA
drop 0.46% and 0.79% in office location, GRL drops 1.08% for the temporal case
due to the intricate structure of minimax game between domain discriminator
and feature extraction. The domain discriminator is unable to entirely
minimize the domain shift while training on fewer data since it does not learn
all feature distinctions of the source and target domains.

CDAN has better results comparing to GAN, GRL, and ADDA, 56.12%, 76.04%, and
92.35% on average of 1-day, 2-day, and 3-day respectively due to its ability
to learn class-wise feature distribution cross-covariance dependency on class
prediction and its feature, therefore, helping understanding the class-
conditional domain drift. However, the accuracy drop happens in 1-day temporal
case which CDAN accuracy is 73.67% while supervised learning method is 76.05%.
The limited amount of training data is still a problem for CDAN as it requires
more data to understand the distribution of source and target.

FixMatch performs best comparing to the other representative methods even
though its original purpose is for semi-supervised learning tasks, not
unsupervised domain adaptation. Because FixMatch synthesizes data utilizing
pseudo-label from target unlabeled data to train the model, it virtually
trains on more data. We observe that the advantage of FixMatch presents more
than the other domain adaptation methods in low data regime as in the temporal
1-day case, FixMatch gives 18.27% performance boost from supervised learning
method while the other methods struggle with adapting for the small dataset.

The average accuracies for 1-day, 2-day, and 3-day cases of GaitSADA are
79.66%, 91.80%, and 97.57% respectively which increase accuracies of vanilla
supervised learning method by a large margin and outperform other baselines.
In particular, our results for the 1-day scenario clearly demonstrate a very
strong accuracy improvement in the low data regime due to the simplicity of
the aligning components without a complex discriminating structure that
requires more training data (e.g., GAN, GRL, ADDA). Specifically, GaitSADA
provides performance boosts by 31.82%, 26.19%, and 17.61% in terms of average
accuracy of 1-day, 2-day, and 3-day cases respectively, compared with the
vanilla supervised learning method, and by 15.41%, 5.83%, and 2.63%, compared
with FixMatch, which is the best baseline.

###  V-D Visualization

Grad-CAM [34] visualizes the weighted feature maps to create a heatmap that
highlights the most important regions to make a prediction. Fig. 6 illustrates
diverse Grad-CAM saliency heatmaps for each algorithm. The other methods’
wrong predictions tend to focus only on a specific part of the spectrogram
data, high-frequency or low-frequency part. The high-frequency part represents
the fast-moving limbs i.e., arms and legs, indicating the importance of this
region for gait recognition. However, (a) Supervised, (b) GAN, and (f)
FixMatch show the models tend to confuse high-frequency and the white noise in
the background. Concentrating on this background noise area and ignoring the
other parts lead to misprediction. (c) GRL, (d) ADDA, and (e) CDAN tend to
focus only on the low-frequency part, the torso region, which may contain less
unique information compared to other parts of the body. We visualize saliency
of GaitSADA on the same data that the other methods predict incorrectly. In
contrast to the other methods, GaitSADA looks broadly at multiple parts of the
body without focusing too much on specific parts which allow our method to
effectively use more unique information from different parts of the body for
gait recognition and thus easily adapt to domain shift.

Figure 6: Grad-CAM [34] visualizations of the wrong prediction of other
methods (a)-(f) comparing to our correct prediction of the same data. The red
color area represents the more important weight contributing to the
prediction.

###  V-E Ablation Study

TABLE III: Ablation experiments of 2-day data on laboratory location as source
data while adapting to different days (i.e., temporal) and 3 different
locations, (i.e., server, conference, and office) using multiple combinations
of Additive Margin, Contrastive learning, Consistency Training and Centroid
Alignment. The results report in % of accuracy.

Stage Components Temporal Server Conference Office Average 1 stage (train
jointly) 1\. supervised 88.46 53.05 71 49.92 65.61 2\. supervised + AM 91.8
73.72 82.8 67.07 78.85 3\. supervised + AM + contrastive learning 94..21 80.38
89.37 70.87 83.71 4\. supervised + AM + consistency + centroid 99.19 52.15
98.2 46.55 74.02 2 stages 5\. (1st stage) supervised + AM + contrastive
learning (2nd stage) supervised + AM + consistency 94.68 78.58 90.3 72.77
84.08 6\. (1st stage) supervised + AM + contrastive learning (2nd stage)
supervised + AM + consistency + centroid 97.63 79.88 99.10 90.59 91.80

We experiment with 2-day data on laboratory location as source data while
adapting to different days (i.e., temporal) and 3 different locations, (i.e.,
server, conference, and office) in the same setting as the experiments in
Table II and examine the behavior of each component contributing to the final
results in both jointly training and 2-stage training setting. We study six
different combinations of components as shown in Table III, i.e., supervised
learning, additive margin softmax loss, contrastive learning, consistency
training, and centroid alignment. In particular, The AM component alone
improves the performance by a large margin. In fact, this component already
outperforms most of the other adversarial-based methods as AM loss increases
inter-class separability and intra-class compactness without training multiple
models (i.e. feature extractor and domain discriminator) as adversarial-based
methods did, which require more data in order to learn the feature
distribution. Contrastive learning further improves by 4.86% on average as the
model learns noise invariance on the target data. Moreover, it can be observed
that two-stage approach outperforms ($84\%$ and $91.8\%$) the single-stage
approach (65.61%, 78.85%, 83.71%, and 74.02%), which indicates the
consistency-based fine-tuning at the second training stage plays an active
role. In addition, it shows that by combining centroid alignment and
consistency training can further improve the recognition accuracy from $84\%$
to $91.8\%$.

##  VI Conclusion

In this paper, we designed a system for gait recognition using a commercial
mmWave radar, conducted multiple times and environments to study the domain
discrepancy in various cases, and proposed a novel domain adaptation method
designed explicitly for mmWave signals. We experiment with five regular domain
adaptation methods to close the domain gaps. The experiments show the
limitations of adversarial-based domain adaptation methods which require more
data in order to learn the data distribution due to the complicated minimax
game of the feature extractor and domain discriminator model. On the other
hand, GaitSADA directly minimizes intra-class variation while maximizing
inter-class deviation by utilizing self- and semi-supervised to learn target
unlabelled data. Our proposed method outperforms standard domain adaptation
methods by a large margin, particularly in the low data regime. Visualization
and ablation studies further validate the effectiveness of GaitSADA.

## References

  * [1] S. Wang, J. Zhang, and L. Tong, “Walk the walk: Attacking gait biometrics by imitation,” in _ISC_ , 2010.
  * [2] C. Fan, Y. Peng, C. Cao, X. Liu, S. Hou, J. Chi, Y. Huang, Q. Li, and Z. He, “Gaitpart: Temporal part-based model for gait recognition,” in _CVPR_ , 2020\.
  * [3] H. Chao, Y. He, J. Zhang, and J. Feng, “Gaitset: Regarding gait as a set for cross-view gait recognition,” _ArXiv_ , vol. abs/1811.06186, 2018.
  * [4] B. Lin, S. Zhang, and X. Yu, “Gait recognition via effective global-local feature representation and local temporal aggregation,” _ICCV_ , pp. 14 628–14 636, 2020.
  * [5] E. Pinyoanuntapong, A. Ali, P. Wang, M. Lee, and C. Chen, “Gaitmixer: Skeleton-based gait representation learning via wide-spectrum multi-axial mixer,” _ArXiv_ , vol. abs/2210.15491, 2022.
  * [6] T. Teepe, A. R. Khan, J. Gilg, F. Herzog, S. Hörmann, and G. Rigoll, “Gaitgraph: Graph convolutional network for skeleton-based gait recognition,” _ICIP_ , pp. 2314–2318, 2021.
  * [7] Y. Zeng, P. H. Pathak, and P. Mohapatra, “Wiwho: Wifi-based person identification in smart spaces,” in _2016 15th ACM/IEEE International Conference on Information Processing in Sensor Networks (IPSN)_ , 2016, pp. 1–12.
  * [8] W. Wang, A. X. Liu, and M. Shahzad, “Gait recognition using wifi signals,” in _UbiComp_ , 2016.
  * [9] A. Pokkunuru, K. Jakkala, A. Bhuyan, P. Wang, and Z. Sun, “Neuralwave: Gait-based user identification through commodity wifi and deep learning,” in _IECON_ , 2018.
  * [10] P. Janakaraj, K. Jakkala, A. Bhuyan, Z. Sun, P. Wang, and M. Lee, “Star: Simultaneous tracking and recognition through millimeter waves and deep learning,” in _2019 12th IFIP Wireless and Mobile Networking Conference (WMNC)_. IEEE, 2019, pp. 211–218.
  * [11] P. Cao, W. Xia, M. Ye, J. Zhang, and J. Zhou, “Radar-id: human identification based on radar micro-doppler signatures using deep convolutional neural networks,” _IET Radar, Sonar & Navigation_, vol. 12, no. 7, pp. 729–734.
  * [12] X. Yang, J. Liu, Y. Chen, X. Guo, and Y. Xie, “Mu-id: Multi-user identification through gaits using millimeter wave radios,” in _IEEE INFOCOM_ , 2020.
  * [13] M. Z. Ozturk, C. Wu, B. Wang, and K. J. Ray Liu, “Gait-based people identification with millimeter-wave radio,” in _IEEE WF-IoT_ , 2021.
  * [14] Q. Xu, Y. Chen, B. Wang, and K. J. R. Liu, “Radio biometrics: Human recognition through a wall,” _IEEE Transactions on Information Forensics and Security_ , vol. 12, no. 5, p. 1141–1155, 2017.
  * [15] F. Adib, C.-Y. Hsu, H. Mao, D. Katabi, and F. Durand, “Capturing the human figure through a wall,” _ACM Transactions on Graphics_ , vol. 34, no. 6, 2015\.
  * [16] J. G. Moreno-Torres, T. Raeder, R. Alaiz-RodríGuez, N. V. Chawla, and F. Herrera, “A unifying view on dataset shift in classification,” _Pattern recognition_ , vol. 45, no. 1, pp. 521–530, 2012.
  * [17] E. Tzeng, J. Hoffman, K. Saenko, and T. Darrell, “Adversarial discriminative domain adaptation,” _CVPR_ , pp. 2962–2971, 2017.
  * [18] Y. Ganin and V. S. Lempitsky, “Unsupervised domain adaptation by backpropagation,” _ArXiv_ , vol. abs/1409.7495, 2015.
  * [19] M. Long, Z. Cao, J. Wang, and M. I. Jordan, “Domain adaptation with randomized multilinear adversarial networks,” _ArXiv_ , vol. abs/1705.10667, 2017.
  * [20] Q. Xie, Z. Dai, E. H. Hovy, M.-T. Luong, and Q. V. Le, “Unsupervised data augmentation for consistency training,” _arXiv: Learning_ , 2019.
  * [21] K. Sohn, D. Berthelot, C.-L. Li, Z. Zhang, N. Carlini, E. D. Cubuk, A. Kurakin, H. Zhang, and C. Raffel, “Fixmatch: Simplifying semi-supervised learning with consistency and confidence,” _ArXiv_ , vol. abs/2001.07685, 2020.
  * [22] D. Berthelot, N. Carlini, I. J. Goodfellow, N. Papernot, A. Oliver, and C. Raffel, “Mixmatch: A holistic approach to semi-supervised learning,” _ArXiv_ , vol. abs/1905.02249, 2019.
  * [23] D. Berthelot, N. Carlini, E. D. Cubuk, A. Kurakin, K. Sohn, H. Zhang, and C. Raffel, “Remixmatch: Semi-supervised learning with distribution alignment and augmentation anchoring,” _ArXiv_ , vol. abs/1911.09785, 2019.
  * [24] T. Chen, S. Kornblith, M. Norouzi, and G. E. Hinton, “A simple framework for contrastive learning of visual representations,” _ArXiv_ , vol. abs/2002.05709, 2020.
  * [25] J.-B. Grill, F. Strub, F. Altch’e, C. Tallec, P. H. Richemond, E. Buchatskaya, C. Doersch, B. Á. Pires, Z. D. Guo, M. G. Azar, B. Piot, K. Kavukcuoglu, R. Munos, and M. Valko, “Bootstrap your own latent: A new approach to self-supervised learning,” _ArXiv_ , vol. abs/2006.07733, 2020.
  * [26] A. van den Oord, Y. Li, and O. Vinyals, “Representation learning with contrastive predictive coding,” _ArXiv_ , vol. abs/1807.03748, 2018.
  * [27] K. He, H. Fan, Y. Wu, S. Xie, and R. B. Girshick, “Momentum contrast for unsupervised visual representation learning,” _CVPR_ , 2019.
  * [28] K. He, X. Zhang, S. Ren, and J. Sun, “Deep residual learning for image recognition,” _CVPR_ , 2015.
  * [29] F. Wang, J. Cheng, W. Liu, and H. Liu, “Additive margin softmax for face verification,” _IEEE Signal Processing Letters_ , vol. 25, pp. 926–930, 2018\.
  * [30] A. Tarvainen and H. Valpola, “Mean teachers are better role models: Weight-averaged consistency targets improve semi-supervised deep learning results,” in _NIPS_ , 2017.
  * [31] S. Laine and T. Aila, “Temporal ensembling for semi-supervised learning,” _ArXiv_ , vol. abs/1610.02242, 2016.
  * [32] I. Goodfellow, J. Pouget-Abadie, M. Mirza, B. Xu, D. Warde-Farley, S. Ozair, A. Courville, and Y. Bengio, “Generative adversarial networks,” _Communications of the ACM_ , vol. 63, no. 11, pp. 139–144, 2020.
  * [33] J. Hoffman, E. Tzeng, T. Darrell, and K. Saenko, “Simultaneous deep transfer across domains and tasks,” _ICCV_ , pp. 4068–4076, 2015.
  * [34] R. R. Selvaraju, A. Das, R. Vedantam, M. Cogswell, D. Parikh, and D. Batra, “Grad-cam: Visual explanations from deep networks via gradient-based localization,” _International Journal of Computer Vision_ , vol. 128, pp. 336–359, 2016.


Pricing the COVID-19 Vaccine: A Mathematical Approach111Preprint submitted on
29 Dec. 2020

Susan Martonosia,222Corresponding author

<EMAIL_ADDRESS>
Banafsheh Behzadb

<EMAIL_ADDRESS>
Kayla Cummingsc

<EMAIL_ADDRESS>
aDepartment of Mathematics, Harvey Mudd College
bDepartment of Information Systems, College of Business Administration at
California State University, Long Beach
cOperations Research Center, Massachusetts Institute of Technology

Abstract

According to the World Health Organization, development of the COVID-19
vaccine is occurring in record time. Administration of the vaccine has started
the same year as the declaration of the COVID-19 pandemic. The United Nations
emphasized the importance of providing COVID-19 vaccines as “a global public
good”, which is accessible and affordable world-wide. Pricing the COVID-19
vaccines is a controversial topic. We use optimization and game theoretic
approaches to model the COVID-19 U.S. vaccine market as a duopoly with two
manufacturers Pfizer-BioNTech and Moderna. The results suggest that even in
the context of very high production and distribution costs, the government can
negotiate prices with the manufacturers to keep public sector prices as low as
possible while meeting demand and ensuring each manufacturer earns a target
profit. Furthermore, these prices are consistent with those currently
predicted in the media.

Keywords: coronavirus, COVID-19, healthcare, optimization, game theory, public
policy, vaccine pricing

##  1 Introduction

The first case of Coronavirus Disease 2019 (COVID-19) was recognized in Wuhan,
China, in December 2019 and was then announced as pandemic by the World Health
Organization (WHO) on March 11th, 2020 ([36]). This pandemic has become one of
the world’s main public health challenges due to both transmission speed and
scale, and the fact that certain segments of the population are at increased
risk of severe illness and death ([30]). As of early December 2020, there have
been more than 66M confirmed cases and 1.5M deaths worldwide due to COVID-19
([41]).

Vaccine development is a complicated process normally requiring years of
research and development. There are several ongoing efforts around the world
to discover a vaccine that mitigates COVID-19 transmission. New-generation
vaccines including messenger RNA (mRNA) vaccines are highly effective and have
potential for fast development, low cost manufacturing and safe
administration. Thus, these vaccines can serve as an assuring substitute for
conventional vaccines ([23]). As of Dec. 5, 2020, the Coronavirus Vaccine
Tracker by New York Times ([20]) reports that there are 58 COVID-19 vaccines
in clinical trials on humans, of which there are two leading manufacturers in
the U.S. in Phase 3 of developing mRNA vaccines: Pfizer-BioNTech and Moderna.
The collaboration between Pfizer and the German company BioNTech reports their
vaccine to be 95% effective in preventing COVID-19 ([26]), and their vaccine
was approved for emergency use in the U.K. on December 2, 2020 ([35]).
(Throughout this work, we refer to the Pfizer-BioNTech collaboration as simply
“Pfizer”). In partnership with National Institute of Health (NIH), Moderna
began vaccine development in January 2020; on Nov. 16, 2020, Moderna announced
its vaccine is 94.5 percent effective, joining Pfizer on the frontier of the
worldwide endeavor to develop a COVID-19 vaccine ([11]).

On July 22, 2020, Pfizer announced an agreement with the U.S. government for
an initial order of 100 million doses of its mRNA vaccine for $1.95 billion
with a possibility of acquiring up to 500 million additional doses ([5]). On
August 11, 2020, Moderna, Inc., announced that the U.S. government committed
up to $2.48 billion for early access to 100 million doses of Moderna’s mRNA
vaccine, with the option to purchase up to 400 million additional doses ([4]).
According to an article published on Nov. 17, 2020, Pfizer has set the initial
price at $19.50 a dose ($39 per patient), and Moderna has set its vaccine
price to $25 per dose ($50 per patient) ([15]). [31] report the price of the
Moderna vaccine as $30 for the two-dose regimen.

Distribution of the mRNA vaccines is a serious challenge due to the need for
ultra-cold storage. The Pfizer and Moderna vaccines require different storage
temperatures: $-70$ degrees Celsius for Pfizer and $-20$ degrees Celsius for
Moderna ([32]). Moderna’s vaccine can be stored in a regular freezer, making
its distribution less costly. Pfizer has provided several ultra cold freezers
to store and distribute their own vaccine ([3]).

The Centers for Disease Control and Prevention (CDC) is the primary public
health organization in the U.S. charged with the responsibility of development
and fulfillment of vaccines and ensuring adequate societal immunization
coverage levels. The Centers for Medicare and Medicaid Services released a
comprehensive plan to ensure all Americans have access to the COVID-19 vaccine
at no cost once available ([6]). According to [31], the COVID-19 vaccine
should be free to everyone in the U.S. under the CDC’s COVID-19 Vaccination
Program ([37]). The CDC discusses vaccine coverage under five program types:
Medicare, Vaccine for Children (VFC), Medicaid, Section 317 (Vaccines for
uninsured adults), the Veterans Administration, and private insurance. We
characterize the U.S. COVID-19 vaccine as a two-sector market, similar to
other vaccine markets. The private sector includes the vaccines that are
purchased by private insurers, and the public sector includes all other
governmental programs. Typically public sector vaccine costs are lower.

Several studies in the operations research literature use theoretical and
empirical approaches to analyze multiple aspects of the COVID-19 pandemic.
These aspects include modeling the spread of the disease ([14]), ICU capacity
management ([1]), importance of wearing high quality masks ([29]), ventilator
inventory management ([16]), the impact of COVID-19 on global supply chain
performance and food supply chain ([13] and [33]), and plasma donations during
the COVID-19 pandemic ([19]), among others. To the best of our knowledge,
there are no theoretical or empirical studies on COVID-19 vaccine pricing.
This study fills this gap.

Previous work by [7] models the negotiation of public sector pediatric vaccine
prices by the CDC in the context of a duopoly. The authors develop an
optimization model for minimizing public sector expenditures while ensuring
adequate public sector vaccine supply and adequate profit margins for each
manufacturer (preventing a monopoly). They treat the private sector market as
an exogenous duopoly price competition operating under the surplus conditions
of a Bertrand-Edgeworth-Chamberlin (BEC) equilibrium, as derived in [2].

This paper considers the U.S. COVID-19 vaccine market as a duopoly between
Pfizer and Moderna. We adapt the model of [7] to determine prices for the U.S.
public sector COVID-19 vaccines. The analysis suggests that even in the
context of high production and distribution costs, the U.S. government can
negotiate prices with the manufacturers that keep public sector prices as low
as possible while meeting demand and ensuring each manufacturer earns a target
profit. Moreover, these prices are consistent with those currently predicted
in the media. Prices are likely to decrease over time as the cold chain
becomes more efficient. Finally, the model presented here can be adapted to
the context of an oligopoly when new vaccines enter the market.

Section 2 presents our adaptation of the model to the COVID-19 vaccine
scenario. In Section 3 we justify the parameter values used for our empirical
tests. Section 4 presents the analysis and discussion. Section 5 concludes.

##  2 Model

We use the same notation and largely the same model as in [7]. The sets,
variables and parameters used in the optimization model can be found in Table
1. The full optimization model is given in Figure 1, which we describe now.

Name | Type |  Meaning |  Estimated Value(s)
---|---|---|---
$M$ | set |  manufacturers |
$S$ | set |  sectors {pub, priv} |
$p_{i}^{\text{pub}}$ | variable |  public sector vaccine price, $i\in M$ |
$q_{i}^{\text{pub}}$ | variable |  public sector vaccine quantity, $i\in M$ |
$z$ | variable |  manufacturers’ absolute public sector price difference |
$p_{i}^{\text{priv}}$ | parameter |  private sector vaccine price, $i\in M$ |  Eq. (9)
$q_{i}^{\text{priv}}$ | parameter |  private sector vaccine quantity, $i\in M$ |  Eq. (10)
$D$ | parameter |  total U.S. demand |  [157.05M, 174.5 M, 191.95M]
$a^{s},\,b,\,c$ | parameters |  demand curve coefficients, $s\in S$ |  Eqs. (14), (11), (12)
$k$ | parameter |  order of magnitude difference between prices and quantities per manufacturer and sector† |  $6$
$U$ | parameter |  private sector capacity lower bound for surplus |  Eq. (13)
$K_{i}$ | parameter |  total production capacity, $i\in M$ |  $K_{Pf}=250$M, $K_{Mod}=200$M
$P_{i}$ | parameter |  minimum profit threshold, $i\in M$ |  $P_{Pf}:\left[\$25.7M,\$234M,\$2.57B\right]$ $P_{Mod}:\left[\$41.4M,\$496M,\$2.57B\right]$
$d_{i}$ | parameter |  two-dose production and distribution cost, $i\in M$ † |  $\left[\$0,\$6.6,\$23.44,\$31.96\right]$
$\gamma$ | parameter |  product similarity |  [0.25, 0.50, 0.75]
$r^{s}$ | parameter |  % total demand sold in sector $s\in S$ † |  $r^{\text{pub}}=0.57$
$\mu$ | parameter |  objective function weight‡ |  0.9
Table 1: Components of the optimization model in Figure 1, and their estimated
base value and range for sensitivity analysis. Quantities, prices and costs
reflect two-dose regimens of the vaccine.
†: new parameters not used in [7]1.
‡: For numerical stability, all vaccine quantities are scaled to units of one
million in the optimization model.
1[7] include a parameter $\rho$ representing a legislative limit on year-over-
year price increases. No such legislation has yet been introduced for the
COVID-19 vaccine ([31]), so we omit this requirement.

$\text{minimize }\quad\mu\Big{(}\sum_{i\in M}q_{i}^{\text{pub}}\cdot
p_{i}^{\text{pub}}\Big{)}+(1-\mu)z$ (1)

$\displaystyle\sum_{i\in M}q_{i}^{\text{pub}}$ $\displaystyle\geq
r^{\text{pub}}D$ (2) $\displaystyle q^{s}_{i}+bp^{s}_{i}-cp^{s}_{j}$
$\displaystyle=a^{s}$ $\displaystyle i,j\in M,\,i\neq j,s\in S$ (3)
$\displaystyle\sum_{s\in S}q^{s}_{i}(p^{s}_{i}-d^{s}_{i})$ $\displaystyle\geq
P_{i}$ $\displaystyle i\in M$ (4) $\displaystyle\sum_{s\in S}q_{i}^{s}$
$\displaystyle\leq K_{i}$ $\displaystyle i\in M$ (5) $\displaystyle
K_{i}-q_{i}^{\text{pub}}$ $\displaystyle\geq U$ $\displaystyle i\in M$ (6)
$\displaystyle\qquad\qquad\qquad z$ $\displaystyle\geq
p_{i}^{\text{pub}}-p_{j}^{\text{pub}}$ $\displaystyle i,j\in M,\,i\neq j$ (7)
$\displaystyle q_{i}^{\text{pub}},\,p_{i}^{\text{pub}},\,z$ $\displaystyle\geq
0$ $\displaystyle i\in M$ (8)

Figure 1:  Model to optimize public-sector price-setting dynamics in a surplus
scenario.

We consider two manufacturers $(i\in M)$ of competing and differentiated
vaccines. The manufacturers sell these vaccines in two sectors, $s\in
S=\\{\text{pub},\text{priv}\\}$. When selling in the private sector the
prices, $p^{\text{priv}}_{i}$, are determined in a Bertrand-Edgeworth-
Chamberlin (BEC) duopoly game presented in [2]. We also assume a surplus
regime, in which each manufacturer’s production capacity, $K_{i}$, for the new
vaccine is sufficiently large that each can maximize private-sector profits
unconstrained by capacity.

Within the public sector, the CDC negotiates quantities, $q^{\text{pub}}_{i}$,
and prices, $p^{\text{pub}}_{i}$, with each manufacturer to minimize the total
government cost while ensuring public sector demand can be fully met. To
reflect the fact that competing vaccines tend to have similar prices, the CDC
chooses among cost-minimizing solutions those which also minimize the absolute
public sector price difference between the two manufacturers,
$z=|p^{\text{pub}}_{1}-p^{\text{pub}}_{2}|$, yielding the weighted objective
function in Equation (1).

The negotiated prices and quantities are constrained by the following
considerations:

  * •

Public sector demand must be met (Eq. 2). We introduce parameter $r^{s}$,
which is the percentage of total demand $D$ attributed to sector $s$, such
that $r^{\text{pub}}+r^{\text{priv}}=1$.

  * •

Prices and quantities follow linear demand curves in each sector (Eq. 3).

  * •

Each manufacturer $i\in M$ must achieve a target profit $P_{i}$ (Eq. 4). As a
novel departure from the model of [7], we introduce per-unit production and
distribution costs $d_{i}$, reflecting the significant expense of the cold
chain required to distribute mRNA vaccines.

  * •

Each manufacturer’s total quantity sold cannot exceed their total capacity
(Eq. 5).

  * •

Each manufacturer’s residual private sector capacity is determined by public
sector quantities and must exceed the bound $U$ governing the BEC game’s
surplus equilibrium regime (Eq. 6).

  * •

The absolute price difference, $z$, between the two manufacturers is expressed
as linear inequality bounds that the objective function drives to tighten (Eq.
7).

The private sector prices and quantities are determined exogenously by the
equilibrium solution to the BEC game given in [7]:

$p_{i}^{\text{priv}}=\frac{a^{\text{priv}}}{2b-c},\forall i\in M.$ (9)
$q_{i}^{\text{priv}}=bp_{i}^{\text{priv}},\forall i\in M,$ (10)

where

$b=\frac{10^{k}}{(1+\gamma)(1-\gamma)},$ (11)

(for $k$ the approximate order of magnitude difference between prices and
quantities produced by each manufacturer in each sector),

$c=\gamma b,$ (12)

and $a^{\text{priv}}$ is estimated as described in Section 3.1. The capacity
bound $U$ is given by:

$U=\frac{a^{\text{priv}}(1+\gamma)}{\gamma}\cdot\left(1-\frac{2(1-\gamma)^{0.5}}{(1+\gamma)^{0.5}(2-\gamma)}\right).$
(13)

Note that Pfizer and Moderna are both anticipated to produce vaccines
requiring a two-dose regimen. We assume all individuals will receive both
doses, so we interpret quantities, unit prices, and unit costs in terms of
individuals served in two-dose increments.

##  3 Parameter Estimation

Here, we discuss estimates and sensitivity analysis ranges for the model
parameters describing the Pfizer and Moderna vaccines, summarized in the last
column of Table 1.

###  3.1 Demand Curve Parameters

We estimate the linear demand curve parameters
$a^{\text{pub}},a^{\text{priv}},b,$ and $c$ in the same manner as [7]. This
method uses publicly available demand and price data to determine the average
public and private sector intercepts, $a^{\text{pub}}$ and $a^{\text{priv}}$.
Historical demand and price data are not yet available for the COVID-19
vaccine, so we use data available for the influenza (flu) vaccine as a proxy.
The flu vaccine is the most widely distributed vaccine in the U.S., with
approximately 175 million U.S. individuals receiving the vaccine in the
2019-20 contract year ([38]). This is the same order of magnitude as the
number of U.S. individuals who could receive the COVID-19 vaccine during its
first year of production. In the U.S. flu vaccine market two primary
manufacturers market vaccines to people older than six months of age: Sanofi
Pasteur (Fluzone®Quadrivalent) and GlaxoSmithKline (Fluarix®Quadrivalent).
Table 2 summarizes public and private sector prices for these two vaccines.

The new COVID-19 vaccines will require two doses per person, as compared to
the flu vaccine’s single dose per person. Additionally, the COVID-19 vaccines
are expected to have much higher distribution costs relative to the flu
vaccine. Thus to estimate $a^{s}$ in each sector, we use the reported number
of administered doses of flu vaccine from the years 2010-11 through 2019-20 to
reflect the number of individuals who will seek the COVID-19 vaccine. We
double the prices given in Table 2 to reflect the cost of a two-dose regimen,
and we add the estimated two-dose production and distribution costs $d_{i}$
that will naturally raise prices. Because two-dose quantities in each sector
by each manufacturer are likely to be on the order of tens of millions, and
prices are likely to be on the order of tens, we set $k=6$. Incorporating
these modifications into the approach of [7], we obtain the following
expression for the estimated linear demand curve intercept:

$a^{s}=\frac{1}{Y}\sum_{y=1}^{Y}\left(\frac{1}{2}\sum_{i\in
M}q_{iy}^{s}+\frac{10^{6}}{2+2\gamma}\sum_{i\in
M}(2p_{iy}^{s}+d_{iy}^{s})\right)\qquad s\in S,$ (14)

where $q_{iy}^{s}$ and $p_{iy}^{s}$ are, respectively, the quantity of flu
vaccine sold and the price set by manufacturer $i\in M$ in sector $s\in S$
during year $y\in\\{1,\dots,Y\\}$. We can interpret $a^{s}$ to be the number
of sector $s$ individuals who would seek a two-dose regimen of the COVID-19
vaccine from either manufacturer if the vaccine had zero cost.

| Public Sector | Private Sector |
---|---|---|---
Year | Sanofi Pasteur | GlaxoSmithKline | Sanofi Pasteur | GlaxoSmithKline | Total Demand (M)
2010-11 | $10.64 | $8.90 | $13.16 | $10.98 | 155.1
2011-12 | $11.68 | $10.97 | $13.16 | $12.41 | 132.0
2012-13 | $11.68 | $9.25 | $14.56 | $10.98 | 134.9
2013-14 | $17.05 | $13.65 | $20.50 | $15.90 | 134.5
2014-15 | $17.44 | $13.65 | $21.09 | $15.90 | 147.8
2015-16 | $17.94 | $14.05 | $21.70 | $16.05 | 146.4
2016-17 | $19.14 | $14.43 | $23.17 | $16.82 | 145.9
2017-18 | $15.68 | $14.43 | $18.72 | $16.82 | 155.3
2018-19 | $15.11 | $13.50 | $19.26 | $16.82 | 169.1
2019-20 | $13.76 | $13.50 | $18.31 | $16.82 | 174.5
Table 2: Yearly flu vaccine prices per dose, by sector ([39]), and total
demand in millions of doses ([38]).

We use the 2019-20 total flu vaccine demand of 174.5M as the baseline two-dose
demand, $D$, considered by our model, and consider values 10% below and above
that baseline.

###  3.2 Production Parameters

The two COVID-19 vaccines are still in development at the time of this
writing, so we rely on information contained in press releases and other media
to estimate production-specific parameters.

####  3.2.1 Total Production Capacity, $K_{i}$.

As stated in Section 1, Pfizer and Moderna have the ability to produce up to
500 million and 400 million doses, respectively, for the U.S. market in 2021.
We use $K_{Pf}=250M$ and $K_{Mod}=200M$ as the capacity measured in two-dose
regimens. Moreover, Pfizer and Moderna have been rapidly ramping up their
supply chains, and it is expected that both manufacturers will have
considerably more capacity than this in 2021, permitting us to safely assume
the surplus regime of the BEC equilibrium. For instance, Pfizer plans to
produce 1.3B doses in 2021 for global use ([3]).

####  3.2.2 Minimum Profit Threshold, $P_{i}$.

Articles in the media question whether Pfizer and Moderna are trying to attain
a profit with the new COVID-19 vaccine ([28, 21]). The U.S. government has an
incentive for the two manufacturers to make a profit on this vaccine, as this
stabilizes the market and also ensures that pharmaceutical companies will be
inclined to once again rapidly develop a vaccine if there is another global
pandemic in the future.

We develop three estimates for each manufacturer’s target profit. First, we
divide each company’s reported annual research and development (R&D) budget by
its reported number of products as an estimate of the average profit it hopes
to recuperate from sales of the new vaccine. Pfizer reports an annual R&D
budget of $8.65B and 337 total pharmaceutical products, corresponding to a
target profit of $25.7M ([24, 25]). However, many of these products are over-
the-counter products or generic labels for other products they carry. If we
instead consider only the 37 products Pfizer lists as its primary products
([27]), we obtain a target profit of $234M. Moderna reports an annual R&D
budget of $496.3M and 12 pharmaceutical products that have reached the
clinical studies phase (none are yet approved), corresponding to a target
profit of $41.4M ([40, 17]). Because Moderna has not yet brought a product to
market, its COVID-19 vaccine could reasonably be expected to recuperate the
entire $496.3M annual R&D budget from the previous year, so we use this value
as a second estimate of Moderna’s target profit. Finally, we can use prior
year revenue estimates using the flu vaccine data reported in Table 2 to
obtain a third estimate of the companies’ target profit, following [7][S.M.
Eq. 7]. In this way, we estimate a target profit of $2.57B for each
manufacturer.

####  3.2.3 Two-Dose Production and Distribution Cost, $d_{i}$.

Per-unit costs associated with the production and distribution of either
manufacturer are not widely known. We use estimates of per-dose investment
costs for the flu vaccine provided by the World Health Organization
([10])[Fig. 2]. In their figure, they report a mean per-dose cost of $\$3.30$
averaged over ten manufacturers and weighted by the quantities each produces.
We can also calculate an unweighted average of $\$15.98$ per dose, and an
unweighted trimmed mean (discarding the smallest and largest reported
investment costs) of $\$11.72.$ We again double these estimates, reflecting a
two-dose regimen of the COVID-19 vaccine, to attain estimates for $d_{i}$ of
$\$0$, $\$6.60$, $\$23.44$, and $\$31.96.$

####  3.2.4 Other Parameters.

The remaining parameters are the product similarity, $\gamma$, the public
sector demand percentage, $r^{\text{pub}}$, and the objective function weight,
$\mu$. The degree of similarity between the two vaccine products depends on
factors such as the mode of delivery and adverse effects. In the absence of
specific information about the vaccines, we test $\gamma=[0.25,0.5,0.75].$ We
assume that the public sector represents $r^{\text{pub}}=57\%$ of the total
U.S. vaccine market ([22]). Lastly, with quantities scaled to millions,
$\mu=0.9$ effectively identifies, among cost-minimizing solutions, a solution
for which the difference in public sector prices is minimal.

##  4 Analysis and Discussion

We tested the 1,296 combinations of model parameters specified in Table 1
using the BARON solver on the NEOS server ([34, 8, 9, 12]). The model is
infeasible whenever both manufacturers have the highest target profit of
$2.57B, regardless of other parameters. The vast majority of remaining
infeasible cases occur when one manufacturer’s target profit is at this
highest value. The high prices required to support such a high target profit
erode demand to the point of failing to meet public sector demand. Of the 554
feasible scenarios, 214 yield public sector prices and quantities that are
nonzero for both manufacturers, reflecting the current expectation that both
manufacturers will participate in the public sector market. We discuss these
214 scenarios.

###  4.1 Prices

Figure 2 plots pairs of public sector prices of Moderna and Pfizer determined
by the model for each of the 214 scenarios, with shading determined by the
total demand used in the scenario.

First, we notice that the prices output by the model are reasonable and
consistent with projections given in the news. The overall centroid of the
data points in Figure 2 is located at a price of $34.39 for Pfizer and $35.09
for Moderna for a two-dose regimen.

Second, we detect three general bands of points, corresponding to the three
values of total demand. When total demand is projected to be low (157.05M two-
dose regimens, in black), optimized prices are generally higher than when
total demand is projected to be high (191.95M two-dose regimens, in light
grey).

Lastly, we notice a general inverse relationship in the prices of the two
manufacturers: higher prices for Pfizer are generally coupled with lower
prices for Moderna, and vice versa. Extreme prices higher than $75 per two-
dose regimen occur when, e.g., both manufacturers have moderately high or high
production costs and product similarity is low, meaning that the manufacturers
have effectively segmented the market to their advantage.

Figure 2: Public sector prices (USD) determined by the optimization model for
two-dose regimens for each manufacturer, as a function of total two-dose
demand (M).

Figure 3 plots the public sector price distribution of each manufacturer by
that manufacturer’s production and distribution cost, $d_{i}$. We see that
public sector prices generally rise with increasing production costs. This is
relevant because cold storage is quite costly ([18]). However, Moderna
recently announced that its vaccine can tolerate a less stringent cold chain
than initially thought, and it is hypothesized that Pfizer might eventually
come to a similar conclusion. Thus, we can expect public sector prices to
start high but eventually decrease as the temperature requirements relax,
lowering distribution costs.

(a) Pfizer public sector price distribution.

(b) Moderna public sector price distribution.

Figure 3: Public sector price distributions determined by the optimization
model for each manufacturer, by each manufacturer’s production and
distribution cost $d_{i}$.

###  4.2 Profit

The CDC has an incentive to ensure minimum target profits are attained by both
manufacturers. This can prevent the formation of a monopoly caused by one
manufacturer leaving the market. Our model identifies public sector prices and
quantities that guarantee each manufacturer achieves its target profit.

News articles have suggested possible prices of the upcoming vaccines by
comparing the federal research funding offered to each manufacturer in
exchange for a contracted number of doses if the vaccine is approved. For
Pfizer, the federal support of $1.95B in exchange for 100M doses indicates an
effective public sector price of $39 per two-dose regimen ([5]). For Moderna,
the federal support of $2.48B in exchange for 100M doses amounts to a public
sector price of $49.6 per two-dose regimen ([4]).

Using our model, we can identify scenarios in which both manufacturers secure
a profit with public sector prices in the ranges cited. Table 3 gives the
settings of the tunable parameters from Table 1 yielding public sector prices
in the range of $\$34-44$ per two doses for Pfizer and $\$45-55$ per two doses
for Moderna. Not only do the realized profits meet the projected target
profits for each manufacturer in each of the five scenarios identified,
realized profit for one of the two manufacturers substantially exceeds target
profit in scenarios 1, 3, and 4.

Scenario | 1 | 2 | 3 | 4 | 5
---|---|---|---|---|---
Total Demand, $D$ (M) | 174.5 | 157.05 | 174.5 | 174.5 | 157.05
Product Similarity, $\gamma$ | 0.75 | 0.75 | 0.75 | 0.75 | 0.75
Pfizer Target Profit, $P_{Pf}$ ($M) | 234 | 234 | 25.7 | 234 | 25.7
Pfizer Production Cost, $d_{Pf}$ ($) | 31.96 | 31.96 | 23.44 | 23.44 | 31.96
Moderna Target Profit, $P_{Mod}$ ($M) | 41.4 | 41.4 | 41.4 | 41.4 | 496
Moderna Production Cost, $d_{Mod}$ ($) | 31.96 | 31.96 | 31.96 | 31.96 | 23.44
Private Sector Equilibrium Price, $p^{\text{priv}}$ ($) | 24.31 | 24.31 | 23.46 | 23.46 | 23.46
Private Sector Equilibrium Quantity, $q^{\text{priv}}$ (M) | 55.6 | 55.6 | 53.6 | 53.6 | 53.6
Pfizer Public Sector Price, $p^{\text{pub}}_{Pf}$ ($) | 40.67 | 40.35 | 38.65 | 38.65 | 38.00
Pfizer Public Sector Quantity, $q^{\text{pub}}_{Pf}$ (M) | 75.7 | 78.6 | 66.7 | 66.7 | 79.7
Pfizer Profit Realized ($M) | 234.0 | 234.0 | 1015.9 | 1015.9 | 25.7
Moderna Public Sector Price, $p^{\text{pub}}_{Mod}$ ($) | 53.66 | 54.92 | 47.15 | 47.15 | 53.86
Moderna Public Sector Quantity, $q^{\text{pub}}_{Mod}$ (M) | 23.8 | 20.3 | 32.7 | 32.7 | 16.3
Moderna Profit Realized ($M) | 90.2 | 41.4 | 41.4 | 41.4 | 496.0
Table 3: Scenarios in which public sector prices fall in the speculated ranges
of $\$34-44$ per two doses for Pfizer and $\$45-55$ per two doses for Moderna.
Bold font indicates realized profit exceeding target profit.

##  5 Conclusion

In summary, we have adapted the model of [7] to the anticipated COVID-19
vaccine market duopoly between Pfizer-BioNTech and Moderna. This paper extends
the original model by incorporating production and distribution costs,
anticipated to be substantial for mRNA vaccines requiring cold storage. Using
publicly available information about the two companies and related vaccine
markets, we use the model to recommend public sector prices to be negotiated
by the CDC under a variety of parameter combinations. Prices yielded by the
model are realistic, averaging $34-$35 per two-dose regimen. Moreover, several
trends emerge:

  1. 1.

The CDC can negotiate prices with the manufacturers that keep public sector
prices as low as possible while meeting demand and ensuring each manufacturer
earns a target profit, even when production costs are high.

  2. 2.

We can expect the public sector costs to decrease over time as the cold
storage requirements of the vaccines loosen.

  3. 3.

If the vaccines are differentiated through factors such as delivery method and
adverse effects, we can expect disparate prices between the two manufacturers.

  4. 4.

Prices speculated in the media are observed as valid outputs of the model,
indicating that such prices could satisfy the public sector demand, achieve
the manufacturers’ target profits, and minimize government costs.

As better information becomes available about the COVID-19 vaccines, the model
can be fine-tuned. Moreover, while this paper has focused on Pfizer-BioNTech
and Moderna, the two frontrunners in the COVID-19 vaccine race, there are
other manufacturers currently developing COVID-19 vaccines intended for the
U.S. market, including several that, as of this writing, are in clinical
trials ([20]). The original game theoretic model of [2] applies more generally
to oligopolies, and the optimization model can also incorporate more than two
manufacturers. As the COVID-19 vaccine market evolves, the CDC and other
decision-makers can use this model to guide pricing decisions.

## References

  * [1] Andres Alban et al.  “ICU capacity management during the COVID-19 pandemic using a process simulation”  In _Intensive Care Medicine_ 46.8, 2020, pp. 1624–1626
  * [2] Banafsheh Behzad and Sheldon Jacobson  “Asymmetric Bertrand-Edgeworth-Chamberlin Competition with Linear Demand: A Pediatric Vaccine Pricing Model”  In _Service Science_ 8.1, 2016, pp. 71–84
  * [3] Karen Bouffard  “Why experts say Moderna COVID-19 vaccine holds advantage over Pfizer’s candidate” Accessed 3 Dec. 2020, Chicago Tribune. https://www.chicagotribune.com/business/ct-biz-covid-19-moderna-pfizer-vaccine-20201117-aoqa3ucfangklgu5megvycfdcm-story.html, 2020
  * [4] Businesswire  “Moderna Announces Supply Agreement with U.S. Government for Initial 100 Million Doses of mRNA Vaccine Against COVID-19 (mRNA-1273)” Accessed 1 Dec. 2020, https://www.businesswire.com/news/home/20200811005852/en/, 2020
  * [5] Businesswire  “Pfizer and BioNTech Announce an Agreement with U.S. Government for up to 600 Million Doses of mRNA-based Vaccine Candidate Against SARS-CoV-2” Accessed 1 Dec. 2020, https://www.pfizer.com/news/press-release/press-release-detail/pfizer-and-biontech-announce-agreement-us-government-600, 2020
  * [6] Centers for Medicare & Medicaid Services  “Trump Administration Acts to Ensure Coverage of Life-Saving COVID-19 Vaccines and Therapeutics” Accessed 4 Dec. 2020, https://www.cms.gov/newsroom/press-releases/trump-administration-acts-ensure-coverage-life-saving-covid-19-vaccines-therapeutics, 2020
  * [7] Kayla Cummings, Banafsheh Behzad and Susan Martonosi  “Centers for Disease Control and Prevention as a Strategic Agent in the Pediatric Vaccine Market: An Analytical Approach”  In _Manufacturing and Service Operations Management_ Articles in Advance, 2020  DOI: https://doi.org/10.1287/msom.2020.0902
  * [8] Joseph Czyzyk, Michael P. Mesnier and Jorge J. Moré  “The NEOS Server”  In _IEEE Journal on Computational Science and Engineering_ 5.3, 1998, pp. 68–75
  * [9] Elizabeth D. Dolan  “The NEOS Server 4.0 Administrative Guide”, 2001
  * [10] Martin Friede et al.  “WHO initiative to increase global and equitable access to influenza vaccine in the event of a pandemic: Supporting developing country production capacity through technology transfer”  In _Vaccine_ 29S, 2011, pp. A2–A7
  * [11] Denise Grady  “Early Data Show Moderna’s Coronavirus Vaccine Is 94.5% Effective” Accessed 16 Nov. 2020, New York Times. https://www.nytimes.com/2020/11/16/health/Covid-moderna-vaccine.html, 2020
  * [12] William Gropp and Jorge J. Moré  “Optimization Environments and the NEOS Server”  In _Approximation Theory and Optimization_ Cambridge University Press, 1997, pp. 167–182
  * [13] Dmitry Ivanov  “Predicting the impacts of epidemic outbreaks on global supply chains: A simulation-based analysis on the coronavirus outbreak (COVID-19/SARS-CoV-2) case”  In _Transportation Research Part E: Logistics and Transportation Review_ , 2020
  * [14] B. Ivorra, M.R. Ferrández, M. Vela-Pérez and A.M. Ramos  “Mathematical modeling of the spread of the coronavirus disease 2019 (COVID-19) taking into account the undetected infections. The case of China”  In _Communications in nonlinear science and numerical simulation_ 88, 2020
  * [15] Katie Jennings  “How Much Will A COVID-19 Vaccine Cost?” Accessed 1 Dec. 2020, Forbes. https://www.forbes.com/sites/katiejennings/2020/11/17/how-much-will-a-covid-19-vaccine-cost/?sh=58c9a732576d, 2020
  * [16] Sanjay Mehrotra et al.  “A model of supply‐chain decisions for resource sharing with an application to ventilator allocation to combat COVID‐19”  In _Naval Research Logistics_ 67.5, 2020  DOI: https://doi.org/10.1002/nav.21905
  * [17] Moderna Inc.  “Company History and Information: mRNA Science - Moderna” Accessed 3 Dec. 2020 at https://www.modernatx.com/about-us/moderna-facts, 2020
  * [18] David Z. Morris  “Pfizer’s COVID vaccine comes with a chilly complication. But that may change” Accessed 6 Dec. 2020, Fortune. https://fortune.com/2020/11/28/pfizer-covid-vaccine-cold-storage-update/, 2020
  * [19] Anna Nagurney and Pritha Dutta  “A Multiclass, Multiproduct COVID-19 Convalescent Plasma Donor Equilibrium Model”, 2020
  * [20] New York Times  “Coronavirus Vaccine Tracker” Accessed 6 Nov. 2020, https://www.nytimes.com/interactive/2020/science/coronavirus-vaccine-tracker.html?auth=login-email&login=email, 2020
  * [21] Janet Nguyen  “How much could Pfizer make from a COVID-19 vaccine?” Accessed 8 Dec. 2020, Marketplace. https://www.marketplace.org/2020/11/12/how-much-could-pfizer-make-from-a-covid-19-vaccine/, 2020
  * [22] Walter Orenstein, Gordon Douglas, Lance Rodewald and Alan Hinman  “Immunizations in the United States: Success, structure, and stress”  In _Health Affairs_ 24.3, 2005, pp. 599–610
  * [23] Norbert Pardi, Michael Hogan, Frederick Porter and Drew Weissman  “mRNA vaccines-a new era in vaccinology”  In _Nature Reviews Drug Discovery_ 17.110493, 2018, pp. 261–279
  * [24] Pfizer Inc.  “2019 Financial Report” Page 16. Accessed 3 Dec. 2020 at https://investors.pfizer.com/financials/annual-reports/default.aspx, 2019
  * [25] Pfizer Inc.  “Our Products Hold the Cure” Accessed 3 Dec. 2020 at https://www.pfizer.com/products, 2020
  * [26] Pfizer Inc.  “Pfizer and BioNTech Conclude Phase 3 study of COVID-19 Vaccine Candidate, Meeting all Primary Efficacy Endpoints” Accessed 8 Dec. 2020 at https://www.pfizer.com/news/press-release/press-release-detail/pfizer-and-biontech-conclude-phase-3-study-covid-19-vaccine, 2020
  * [27] Pfizer Inc.  “Pfizer Company Fact Sheet” Accessed 3 Dec. 2020 at https://www.pfizer.com/about/leadership-and-structure/company-fact-sheet, 2020
  * [28] Shawn Radcliffe  “How Much Will You Pay for a COVID-19 Vaccine? Here’s What We Know” Accessed 8 Dec. 2020, Healthline. https://www.healthline.com/health-news/how-much-will-you-pay-for-a-covid-19-vaccine-heres-what-we-know, 2020
  * [29] Erik Rosenstrom et al.  “High-Quality Masks Can Reduce Infections and Deaths in the U.S.”  In _medRxiv_ Cold Spring Harbor Laboratory Press, 2020  DOI: 10.1101/2020.09.27.20199737
  * [30] Wagner Gouvea Dos Santos  “Natural history of COVID-19 and current knowledge on treatment therapeutic options”  In _Biomed Pharmacother_ 129.110493, 2020
  * [31] Karyn Schwartz et al.  “Vaccine Coverage, Pricing, and Reimbursement in the U.S.”, 2020
  * [32] Selena Simmons-Duffin  “Why Does Pfizer’s COVID-19 Vaccine Need To Be Kept Colder Than Antarctica?” Accessed 2 Dec. 2020, NPR. https://www.npr.org/sections/health-shots/2020/11/17/935563377/why-does-pfizers-covid-19-vaccine-need-to-be-kept-colder-than-antarctica, 2020
  * [33] Sube Singh, Ramesh Kumar, Rohit Panchal and Manoj Kumar Tiwari  “Impact of COVID-19 on logistics systems and disruptions in food supply chain”  In _International Journal of Production Research_ 0.0 TaylorFrancis, 2020, pp. 1–16  DOI: 10.1080/00207543.2020.1792000
  * [34] M. Tawarmalani and N. V. Sahinidis  “A polyhedral branch-and-cut approach to global optimization”  In _Mathematical Programming_ 103, 2005, pp. 225–249
  * [35] Katie Thomas  “New Pfizer Results: Coronavirus Vaccine Is Safe and 95% Effective” Accessed 8 Dec. 2020, https://www.nytimes.com/2020/11/18/health/pfizer-covid-vaccine.html, 2020
  * [36] United States Centers for Disease Control and Prevention  “About COVID-19” Accessed 6 Nov. 2020, https://www.cdc.gov/coronavirus/2019-ncov/cdcresponse/about-COVID-19.html, 2020
  * [37] United States Centers for Disease Control and Prevention  “COVID-19 Vaccination Program Interim Playbook for Jurisdiction Operations” Accessed 8 Dec. 2020, https://www.cdc.gov/vaccines/imz-managers/downloads/COVID-19-Vaccination-Program-Interim_Playbook.pdf, 2020
  * [38] United States Centers for Disease Control and Prevention  “Historical Reference of Seasonal Influenza Vaccine Doses Distributed” Accessed 19 Nov. 2020 at https://www.cdc.gov/flu/prevent/vaccine-supply-historical.htm., 2020
  * [39] United States Centers for Disease Control and Prevention  “VFC CDC Vaccine Price List Archives” Accessed 24 Nov. 2020 at https://www.cdc.gov/vaccines/programs/vfc/awardees/vaccine-management/price-list/archive.html., 2020
  * [40] Wall Street Journal  “MRNA — Moderna Inc. Annual Income Statement - WSJ” Accessed 3 Dec. 2020 at https://www.wsj.com/market-data/quotes/MRNA/financials/annual/income-statement, 2019
  * [41] World Health Organization  “WHO Coronavirus Disease (COVID-19) Dashboard” Accessed 8 Dec. 2020, https://covid19.who.int/, 2020

Acknowledgement. The third author is supported by the National Science
Foundation Graduate Research Fellowship under Grant No. 1122374. The opinion,
findings, and conclusions or recommendations expressed in this material are
those of the authors and do not necessarily reflect the opinions of the
National Science Foundation.