${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{D}^{\text{User}}}}=\\{\mathbf{D}\in\mathbb{R}^{n\times
k}:0\leq\|\mathbf{D}[i,:]\|_{1}\leq r_{i},1\leq i\leq n\\}$, then the square
of sensitivity of $\mathbf{B}$ is,

$\displaystyle\max_{\mathbf{D}\in{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{D}^{\text{User}}}}}\sum_{1\leq
i\leq n}\frac{\|\mathbf{D}[i,:]\|_{2}^{2}}{\sigma_{i}^{2}}$
$\displaystyle\leq\sum_{1\leq i\leq
n}\frac{\|\mathbf{D}[i,:]\|_{1}^{2}}{\sigma_{i}^{2}}$
$\displaystyle\leq\sum_{1\leq i\leq n}\frac{r_{i}^{2}}{\sigma_{i}^{2}}.$

∎

###  8.5. Proof for Theorem 4.7

###### Proof.

Let $f(\mathbf{X})=\mathbf{B}\mathbf{X}$, from theorem 4.6 the sensitivity of
$f$ is
$\Delta_{2}(f)=\sqrt{\frac{r_{1}^{2}}{\sigma_{1}^{2}}+\cdots+\frac{r_{n}^{2}}{\sigma_{n}^{2}}}$.
The constraint in problem 5 makes sure that $\Delta_{2}^{2}(f)\leq 2\rho$. We
can flatten $f(\mathbf{X})$ and $\mathbf{Z}\sim
MVG_{n,k}(0,\mathbf{I}_{n},\mathbf{I}_{k})$ as vectors of length $n\times k$.
From Theorem 2.5, we get
$f(\mathbf{X})+\mathbf{Z}=\mathbf{B}\mathbf{X}+\mathbf{Z}$ satisfies
$\rho$-zCDP. Using the post-processing property, we know that
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}(\mathbf{X})$
=
$\mathbf{L}(\mathbf{B}\mathbf{X}+\mathbf{Z})=\mathbf{W}\mathbf{X}+\mathbf{L}\mathbf{Z}$
satisfies $\rho$-zCDP. ∎

###  8.6. Truncation Algorithm

Input: The maximum user contributions on a given publisher at the current day
is $r$. The conversions in day $i$ are
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{c}}}_{1},\cdots,{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{c}}}_{k}$,
the corresponding attributed weights to the given publisher are
$w_{1},\cdots,w_{k}$ (they are either 0 or 1). Suppose there is a function
that maps each conversion to its user, denoted as get_user.

Output: $x$: The truncated advertising measurement on the current day.

1 user_contribution = {}

2 $x=0$

3 for  _$t=1\to k$_ do

4 user =
get_user$({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{c}}}_{t})$

5 if  _user_contribution[user] $\geq r$_ then

6 $w_{t}=0$

7

8 $x$ += $w_{t}$

9 user_contribution[user] += $w_{t}$

10

11return $x$

Algorithm 2 Truncate Data Stream with User Contribution Constraints

Algorithm 2 shows how to truncate the data stream with user contribution
constraints.

###  8.7. Details of Baseline Methods

We compare our algorithm (named Ads-BPC) with other baselines. Unless
otherwise specified, we allocate a total privacy budget of $\rho_{total}=1$.

IPA (case2023interoperable, ): The Interoperable Private Attribution (IPA)
protocol leverages a combination of multi-party computation and differential
privacy to facilitate the processing of individuals’ data while ensuring that
only aggregated measurements are disclosed. In this paper, our primary focus
is on comparing how IPA achieves differential privacy, assuming the presence
of a trusted third party that provides the joined table.

IPA adds independent noise to the data vector
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$,
here the function is
$f({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})={{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$.
Under the event level adjacency
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Event}}}}$,
the $L_{2}$ sensitivity of $f$ is

$\displaystyle\Delta_{2}(f)=$
$\displaystyle\max_{({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}},\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}})\in{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Event}}}}}\|f({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})-f(\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}})\|_{2}$
$\displaystyle=$
$\displaystyle\max_{({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}},\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}})\in{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Event}}}}}\|{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}-\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color<EMAIL_ADDRESS>
The final equation holds because
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$
is derived from
$\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}}$
through a record substitution. This permits modifications in two entries of
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$,
with each entry being altered by a maximum of 1. Consequently, the largest
possible change in the $L_{2}$ norm is $\sqrt{2}$. The mechanism is
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})=\mathbf{W}({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}+N(0,\sigma^{2}\mathbf{I}_{n}))$.
From Theorem 2.5,
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}+N(0,\sigma^{2}\mathbf{I}_{n})$
satisfies $\rho$-zCDP when
$\rho=\frac{\Delta^{2}_{2}(f)}{2\sigma^{2}}=\frac{1}{\sigma^{2}}$, so
$\sigma^{2}=\frac{1}{\rho}$. By the post-processing property 2.6,
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}$
satisfies $\rho$-zCDP. Using the group privacy property 2.8, under user-level
neighboring datasets with a global contribution constraint of at most $GS$,
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}$
satisfies $\rho$-zCDP when $\sigma^{2}=\frac{GS^{2}}{\rho}$.

Id-BPC (add independent noise with bounded per-day contribution): We set a
maximum number of contributions $r$ a user can make in 1 day, if any user
contributes more than $r$ records (conversions) in 1 day, only the first $r$
records will be kept in the dataset, other records are cropped. Then we add
independent noise to the data vector
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$.
The mechanism has the form
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})=\mathbf{L}(\frac{1}{\sigma}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}+\mathbf{z})=\mathbf{W}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}+\mathbf{L}\mathbf{z}$.
Here $\mathbf{z}\sim N(0,\mathbf{I}_{n})$, the strategy matrix is a scaled
version of the identity matrix $\frac{1}{\sigma}\mathbf{I}_{n}$,
$\mathbf{L}=\sigma\mathbf{W}$. Let
$f({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})=\frac{1}{\sigma}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$,
under the user-level neighboring dataset with bounded per-day contributions
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{User}}}}$,
the sensitivity of $f$ is given by

$\displaystyle\Delta_{2}(f)=$
$\displaystyle\max_{({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}},\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}})\in{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{User}}}}}\frac{1}{\sigma}\|{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}-\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}}\|_{2}$
$\displaystyle=$
$\displaystyle\max_{\mathbf{d}\in{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{D}^{\text{User}}}}}\frac{1}{\sigma}\|\mathbf{d}\|_{2}=\frac{\sqrt{n}r}{\sigma},$

where
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{D}^{\text{User}}_{ex}}}=\\{\mathbf{d}\in\mathbb{R}^{n}:d_{i}\in\\{-r,r\\},1\leq
i\leq n\\}$, the last maximum is achieved when $d_{i}=r$. To satisfy
$\rho$-zCDP, $\Delta_{2}^{2}(f)=\frac{nr^{2}}{\sigma^{2}}\leq 2\rho$. We get
that the scale parameter is $\sigma^{2}=\frac{nr^{2}}{2\rho}$.

StreamUDP (dong2023continual, ): StreamUDP adds independent noise to a binary
tree structure. Let $\mathbf{C}$ be the matrix that represents the tree
structures in (dong2023continual, ), for example, when $n=7$,
$\mathbf{C}=diag\\_block(\mathbf{C}_{1},\mathbf{C}_{2},\mathbf{C}_{4})$ with
$\mathbf{C}_{2^{h}}$ being the matrix representation for a complete tree with
height $h$ and

$\displaystyle\mathbf{C}_{1}=\begin{bmatrix}1\end{bmatrix},\mathbf{C}_{2}=\begin{bmatrix}1&0\\\
0&1\\\ 1&1\end{bmatrix},\mathbf{C}_{4}=\begin{bmatrix}1&0&0&0\\\ 0&1&0&0\\\
1&1&0&0\\\ 0&0&1&0\\\ 0&0&0&1\\\ 0&0&1&1\\\ 1&1&1&1\end{bmatrix}.$

Let
$f({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})=\mathbf{C}_{2^{h}}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$.
Under the event level adjacency
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Event}}}}$,
the $L_{2}$ sensitivity of $f$ is

$\displaystyle\Delta_{2}(f)=$
$\displaystyle\max_{({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}},\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}})\in{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Event}}}}}\|f({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})-f(\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}})\|_{2}$
$\displaystyle=$
$\displaystyle\max_{({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}},\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}})\in{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Event}}}}}\|\mathbf{C}_{2^{h}}({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}-\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}})\|_{2}=\sqrt{2h}.$

The final equation remains valid because
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$
is generated from
$\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}}$
through a record substitution. This allows for alterations in two entries
within
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$,
where each entry is subject to a maximum change of 1. Modifying a single entry
in
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$
by 1 has the potential to affect the output of $\mathbf{C}_{2^{h}}$ by $h+1$
entries, as each entry can change by 1. When two entries in
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$
are modified, the maximum potential impact on the output of
$\mathbf{C}_{2^{h}}$ is limited to at most $2(h+1)-2$ entries. The $-2$ term
accounts for the fact that both changes will influence the root. Consequently,
the greatest possible change in the $L_{2}$ norm is $\sqrt{2h}$. The mechanism
is
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})=\mathbf{C}_{2^{h}}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}+N(0,\sigma^{2}\mathbf{I}_{n})$.
From Theorem 2.5,
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}$
satisfies $\rho$-zCDP when $\sigma^{2}=\frac{h}{\rho}$. Using the group
privacy property (Theorem 2.8), under user-level neighboring datasets with a
global contribution constraint of at most $GS$,
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}$
satisfies $\rho$-zCDP when $\sigma^{2}=\frac{h*GS^{2}}{\rho}$. Note that in
(dong2023continual, ), the authors also proposed a method for dynamically
setting the truncation threshold. They allocate a portion of the privacy
budget to monitor the current maximum contributions of users. Whenever a user
surpasses the threshold, they double the threshold. We believe that a similar
approach could be adapted for Ads-BPC. Instead of establishing a predefined
user contribution bound $r$, a portion of the privacy budget could be
allocated to dynamically adjust this bound. We consider it an interesting
direction for future research, which stands apart from the methods discussed
in this paper. To ensure a fair comparison, we maintain a fixed threshold (GS)
for StreamUDP.

Bin-BPC (binary tree method with bounded per-day contribution): We set a
maximum number of contributions $r$ a user can make in 1 day, if any user
contributes more than $r$ records (conversions) in 1 day, only the first $r$
records will be kept in the dataset, other records are cropped. Then we add
independent noise to a binary tree structure. Similar to StreamUDP, let
$\mathbf{C}$ be the matrix that represents the tree structures in
(dong2023continual, ), the mechanism has the form
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})=\mathbf{L}(\frac{1}{\sigma}\mathbf{C}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}+\mathbf{z})=\mathbf{W}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}+\mathbf{L}\mathbf{z}$.
Here $\mathbf{z}\sim N(0,\mathbf{I}_{n})$, the strategy matrix is a scaled
version of $\mathbf{C}$: $\frac{1}{\sigma}\mathbf{C}$, $\mathbf{L}$ represents
how to reconstruct $\mathbf{W}$ from $\frac{1}{\sigma}\mathbf{C}$. Let
$f({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})=\frac{1}{\sigma}\mathbf{C}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$,
under the user-level neighboring dataset with bounded per-day contributions
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{User}}}}$,
the sensitivity of $f$ is

$\displaystyle\Delta_{2}(f)=$
$\displaystyle\max_{({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}},\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}})\in{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{User}}}}}\frac{1}{\sigma}\|\mathbf{C}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}-\mathbf{C}\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}}\|_{2}$
$\displaystyle=$
$\displaystyle\max_{\mathbf{d}\in{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{D}^{\text{User}}}}}\frac{1}{\sigma}\|\mathbf{C}\mathbf{d}\|_{2}=\frac{\sqrt{\sum_{i,j}(\mathbf{C}^{T}\mathbf{C})[i,j]}*r}{\sigma},$

in which,
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{D}^{\text{User}}_{ex}}}=\\{\mathbf{d}\in\mathbb{R}^{n}:d_{i}\in\\{-r,r\\},1\leq
i\leq n\\}$, the last maximum is achieved when $d_{i}=r$, because elements in
$\mathbf{C}^{T}\mathbf{C}$ are all positive. To satisfy $\rho$-zCDP,
$\Delta_{2}^{2}(f)=\frac{\sum_{i,j}(\mathbf{C}^{T}\mathbf{C})[i,j]*r^{2}}{\sigma^{2}}\leq
2\rho$. We get that the scale parameter is
$\sigma^{2}=\frac{\sum_{i,j}(\mathbf{C}^{T}\mathbf{C})[i,j]*r^{2}}{2\rho}$.

UMM (user-level matrix mechanism): It’s adapted from the event-level matrix
mechanism (henzinger2022constant, ). The event-level matrix mechanism the form
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{M}}}({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}})=\mathbf{L}(\mathbf{S}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}+\mathbf{z})$,
where $\mathbf{z}\sim N(0,\mathbf{I}_{n})$ and the matrix $\mathbf{S}$ is
derived by solving an optimization problem. The original matrix mechanism
problem (li2013optimal, ; mckenna2021hdmm, ) works for event-level neighboring
datasets where
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$
is derived from
$\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}}$
by adding or removing a single record. This is different from our event-level
neighboring datasets
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Event}}}}$,
where
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}$
is derived from
$\tilde{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathbf{x}}}}$
by substituting a single record. We can get the matrix mechanism under
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Event}}}}$
by proper scaling. Suppose $\mathbf{S}$ is the optimal solution derived by
solvers in (mckenna2021hdmm, ) under $\rho$-zCDP, then $\mathbf{S}/2$ works
under
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Event}}}}$.
Using the group privacy property, we know that $\mathbf{S}/(2GS)$ works under
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{N}^{\text{Global}}}}$
with a global contribution constraint of at most GS.

###  8.8. Extra Experiments

Figure 4. Comparison of different days $n$ on Facebook Sales Conversion
Dataset. Figure 5. Comparison of different privacy budget $\rho$ on Facebook
Sales Conversion Dataset. Figure 6. Comparison of different user contributions
$r$ on Facebook Sales Conversion Dataset.

In this section, we demonstrate how different parameters influence the errors.
We use the Facebook Sales Conversion dataset for illustration.

####  8.8.1. Experiments on Different Ad Campaign Durations

In this section, we study the performance of different algorithms with respect
to the number of days $n$, $n$ is chosen from [7, 15, 31, 63, 127, 255]. The
dataset we use is the Facebook Sales Conversion dataset with LTA as described
in Section 6.4. The privacy budget is set as $\rho=1$. The user contribution
limit for each day is $r=\lceil\frac{60}{n}*2\rceil$. Figure 4 illustrates the
performance of various algorithms. Ads-BPC outperforms others consistently.

####  8.8.2. Experiment on Different Privacy Budgets

In this section, we study the performance of different algorithms with respect
to the privacy budget $\rho$, $\rho$ is chosen from [0.25, 0.5, 1, 2, 4, 8].
The dataset we use is the Facebook dataset with LTA as described in section
6.4. The number of days is $n=31$. The user contribution limit for each day is
$r=\lceil\frac{60}{31}*2\rceil=4$. As depicted in Figure 5, the errors
associated with various algorithms decrease as $\rho$ increases. This observed
trend can be attributed to the fact that with a larger privacy budget, the
amount of noise added to the dataset decreases, leading to lower error rates.

####  8.8.3. Experiments on Different User Contributions

In this section, we study the performance of different algorithms with respect
to the user contribution limits $r$, $r$ is chosen from [1, 2, 3, 4, 5, 6, 8,
10, 12]. The dataset we use is the Facebook dataset with LTA as described in
Section 6.4. The number of days is $n=31$. The privacy budget is $\rho=1$. In
Figure 6, as the user contribution limit increases, more noise is introduced
into the data, resulting in higher variance. Notably, Ads-BPC demonstrates
improved performance when the threshold $r$ is chosen reasonably.


11affiliationtext:  University of Florida, Warrington College of
Business22affiliationtext:  University of Texas at Austin, McCombs School of
Business

# Nonparametric Discrete Choice Experiments with
Machine Learning Guided Adaptive Design

Mingzhang Yin  Ruijiang Gao  Weiran Lin  Steven M. Shugan

###### Abstract

Designing products to meet consumers’ preferences is essential for a
business’s success. We propose Gradient-based Survey (GBS), a discrete choice
experiment for multiattribute product design. The experiment elicits consumer
preferences through a sequence of paired comparisons for partial profiles. GBS
adaptively constructs paired comparison questions based on the respondents’
previous choices. Unlike the traditional random utility maximization paradigm,
GBS is robust to model misspecification by not requiring a parametric utility
model. Cross-pollinating the machine learning and experiment design, GBS is
scalable to products with hundreds of attributes and can design personalized
products for heterogeneous consumers. We demonstrate the advantage of GBS in
accuracy and sample efficiency compared to the existing parametric and
nonparametric methods in simulations.

##  1 Introduction

Identifying an optimal product based on consumers’ preferences is essential
for the success of a business. Such a problem is prevalent in companies where
the product consists of multiple attributes such as health insurance, cell
phone plans, pizzas, automobiles, logos, and email advertisements
(Balakrishnan et al., , 2004; Netzer and Srinivasan, , 2011; Bertsimas and
Mišić, , 2017; Ellickson et al., , 2022). In this paper, we focus on products
with discrete attributes represented as multivariate binary variables, as is
the case in A/B testing 111A discrete attribute with more than two levels can
be coded as multiple binary variables..

Selecting appropriate products from a choice set is complicated for decision-
makers. First, the desired product should meet the unobserved consumer
preferences. The latent preferences are often revealed by survey experiments
known as conjoint analysis (Luce and Tukey, , 1964; Green and Rao, , 1971).
Since the milestone work of Louviere and Woodworth,  (1983), the choice-based
conjoint (CBC) analysis become one of the most widely used methods to quantify
multiattribute preference (Hein et al., , 2020). A prevalent assumption is
that the preference for an attribute is quantified by a part-worth score, and
the total utility of a product profile is the sum of part-worths (Green and
Srinivasan, , 1978). This parametric assumption, however, may oversimplify how
respondents encode and evaluate products (Allenby et al., , 2005). Second, the
choice set for product design grows exponentially with the number of
attributes, making the problem NP-hard (Kohli and Krishnamurti, , 1989). The
problem is more evident with the development of technology when more and more
components are integrated into a single product. For example, the design of a
smartphone might need to consider hundreds of attributes from a digital
camera, screen display, connectivity modules, software applications, and a
range of sensors. The high-dimensional attributes pose a scalability challenge
to the extant product design methods. Lastly, consumer preference is
heterogeneous, and it is desired to provide a customized product aligned with
individual tastes. Nowadays, an increasing number of products are presented as
digital content, making personalized product design feasible. For example, the
email campaign content can be designed based on a receiver’s demographic
information. The fast development of generative AI might expand the need for
personalized product design.

Figure 1: Demonstration of a paired choice question for a logo design.

We propose GBS that is robust to model misspecification of the utility
function, scalable to high-dimensional attributes and applies to single or
personalized product design. The idea of GBS is to combine gradient-based
machine learning with discrete choice experiments (DCEs). We model the product
attributes as random variables following Bernoulli distributions, then
maximize an objective function such as the market share. The inference of the
optimal product is conducted by computing the gradient of the objective
function with respect to the Bernoulli distribution parameters. However, it is
challenging to compute the gradient because the genuine functional form of
respondent choice given an item is not explicit and the gradient cannot be
directly computed with non-continuous variables. To address these challenges,
we adopt the score function method, a.k.a. the policy gradient method in
reinforcement learning, to compute the gradient using data from DCEs. We build
on the recent development of discrete optimization in machine learning to
develop new variance reduction tools for the score function gradient. The
unbiased and low-variance gradient is then used to generate the next survey
questionnaire and update the optimal product’s distribution parameter by
stochastic gradient descent (SGD).

Bridging gradient-based machine learning and experiments borrows strength from
both worlds. The sequential nature of SGD naturally provides an adaptive
approach to design experiments. Unlike traditional heuristic or rule-based
adaptive design (Green et al., , 1991; Netzer and Srinivasan, , 2011), the
gradient method is derived mathematically and is readily to incorporate new
statistical tools for variance control. The proposed GBS maximizes the
information extracted from each paired comparison question by utilizing the
greedy property and the variance reduction technique of the gradient. GBS
collects the data in the same way as a classic paired conjoint design where
the respondents are randomly selected from a population and asked to choose
between a pair of items with different attributes. Hence, GBS can be
seamlessly integrated into commercial adaptive conjoint software like Sawtooth
Software (Huber, , 2005). On the other hand, the experiments generate the data
in a way by which machine learning can explore combinatorial actions. In
contrast, using observational data may face the overlapping and extrapolation
problems because the data in practice may only be collected on a small subset
of attribute combinations from the large action space.

In simulation, we demonstrate that GBS scales to hundreds of attributes
efficiently and can infer the optimal products accurately. We compare GBS with
parametric and nonparametric baseline methods. Over a variety of utility
functions, GBS is more robust in model specification than parametric methods
and is more sample-efficient than neural networks. Finally, we apply GBS to
learn an individual policy with combinatorial actions via experiments.

##  2 Problem Set Up

Denote $Y_{i}(Z,Z_{0})$ as the user $i$’s choice of product given the focal
product $Z$ and a baseline product $Z_{0}$. The products are represented by
$K$ binary features $Z,Z_{0}\in\\{0,1\\}^{K}$. The baseline product $Z_{0}$
can be a competitor’s product, the exisitng product in the market, or an empty
set. The potential outcome $Y_{i}(Z,Z_{0})$ is distributed according to the
unknown data generating process. $Y_{i}(Z_{1},Z_{0})=1$ if user $i$ chooses
$Z_{1}$ and $Y_{i}(Z_{1},Z_{0})=0$ if $Z_{0}$ is chosen.

We first consider how to identify a single optimal product. The objective is
which represents the potential market share of the proposed product $Z$ in the
presence of the baseline product. The difficulties in solving LABEL:eq:single
are that the number of possible products grows exponentially with feature
numbers, and the functional form of $Y_{i}(Z,Z_{0})$ is unknown and might be
complicated.

##  3 Gradient-based Survey Design

We adopt the Random Utility Maximisation (RUM) framework to model the choice
behavior (McFadden, , 1974). Each product $Z_{j}$ is associated with a utility
$U_{i}(Z_{j})$ for individual $i$. The alternative $Z_{1}$ is chosen from the
pair $(Z_{1},Z_{0})$ if and only if $U_{i}(Z_{1})>U_{i}(Z_{0})$. In RUM, the
utility is decomposed as $U_{i}(Z_{j})=V_{i}(Z_{j})+\epsilon_{ij},j=0,1$ where
$V_{i}(Z_{j})$ is often called the representative utility (Train, , 2009). We
assume the random $\epsilon_{ij}$ independently follows type I extreme value
(Gumbel) distribution. Accordingly, the probability of choosing item $Z$ is
$p(Y_{i}(Z,Z_{0})=1)=\exp(V_{i}(Z))/(\exp(V_{i}(Z))+\exp(V_{i}(Z_{0})))$.
Notice that the product $Z$ that maximizes the choice probability
$p(Y_{i}(Z,Z_{0})=1)$ is irrelevant to what the baseline product $Z_{0}$ is.
We do not make a parametric assumption for the representative utility
$V_{i}(Z)$.

Instead of optimizing the discrete features $Z$ directly, we consider $Z$ as
random variables following distribution
$p(Z;\pi)=\prod_{k=1}^{K}\text{Bern}(z_{k};\pi_{k})$, $\pi_{k}\in[0,1]$. We
then transform the problem in LABEL:eq:single to an equivalent problem with
the same optimal solution as

$\displaystyle\max_{\pi}V(\pi)=\mathbb{E}_{Z\sim
p(Z;\pi)}\mathbb{E}[Y_{i}(Z,Z_{0})\,|\,Z].$ (1)

The equivalence of probabilistic reformulation is shown in Yin et al.,  (2020)
Theorem 1. To facilitate the optimization in an unconstrained space, we
parameterize the probability by the sigmoid function
$Z\sim\prod_{k=1}^{K}\text{Bern}(z;\pi_{k}=\sigma(\phi_{k})),\phi\in\mathbb{R}^{K}$,
$\sigma(x)=1/(1+e^{-x})$, and optimize $V(\phi):=V(\pi=\sigma(\phi))$ with the
logits $\phi$.

The gradient of Eq. 1 can be computed by the score function estimator (a.k.a.
REINFORCE) as $\nabla_{\phi}V(\phi)=\mathbb{E}_{Z\sim
p(Z;\phi)}[\nabla_{\theta}\log p(Z;\phi)\mathbb{E}[Y_{i}(Z,Z_{0})\,|\,Z]]$.
However, score function gradients often suffer from high variance, and many
works have been devoted to reducing the variance. With a direct application of
the antithetic sampling and control variates (Yin and Zhou, , 2019) to the
product design problem, we have the following Lemma.

######  Lemma 1.

An unbiased gradient of the objective $V(\phi)$ is

$\displaystyle\nabla_{\phi}V(\phi)=\mathbb{E}_{u\sim\prod_{k=1}^{K}\text{Unif}(0,1)}\Big{[}\mathbb{E}\big{[}(Y_{i}(Z_{1}(u),Z_{0})-Y_{i}(Z_{2}(u),Z_{0}))(u-\frac{1}{2})\,|\,u\big{]}\Big{]},$
(2)

where
$Z_{1}(u)=\mathbf{1}[u>\sigma(-\phi)],Z_{2}(u)=\mathbf{1}[u<\sigma(\phi)].$

For completeness, the proof is in the Appendix. The gradient in Eq. 2 is
guaranteed to reduce the variance of score function gradient for non-negative
objective (Yin and Zhou, , 2019). A Monte-Carlo estimate of the gradient in
Eq. 2 can be computed by asking respondent $i$ to choose between products
$Z_{1}$ and $Z_{0}$, between $Z_{2}$ and $Z_{0}$, and take the difference of
the choices. However, such question design might be sensitive to the selection
of the baseline $Z_{0}$. A strong or weak $Z_{0}$ may make the choices
$Y_{i}(Z_{1},Z_{0})$ and $Y_{i}(Z_{2},Z_{0})$ the same frequently, resuting in
a zero gradient and slow learning. Moreover, a respondent’s preferences across
two choices may not be fully consistent and comparable.

To mitigate these problems, we marginalize out the zero gradients using Lemma
2.

######  Lemma 2.

$Y(Z_{1},Z_{0})-Y(Z_{2},Z_{0})\,|\,Y(Z_{1},Z_{0})\neq
Y(Z_{2},Z_{0})\stackrel{{\scriptstyle d}}{{=}}(2Y(Z_{1},Z_{2})-1)$

With Lemmas 2 and 1, we have the following result.

######  Lemma 3.

$\tilde{g}=(2Y_{i}(Z_{1}(u),Z_{2}(u))-1)(u-\frac{1}{2})p(\mathcal{A}_{i})$
satisfies $\mathbb{E}[\tilde{g}]=\nabla_{\phi}V(\phi)$, where $Z_{1}(u)$,
$Z_{2}(u)$ are defined as in Lemma 1,$u\sim\prod_{k=1}^{K}\text{Unif}(0,1)$,
event $\mathcal{A}_{i}=\\{Y_{i}(Z_{1},Z_{0})\neq Y_{i}(Z_{2},Z_{0})\\}$.

Computing $\tilde{g}$ takes a choice $Y_{i}(Z_{1}(u),Z_{2}(u))$ from a random
respondent $i$ between two partial profiles $Z_{1}(u)$ and $Z_{2}(u)$. The
baseline product $Z_{0}$ is no longer needed. The partial profile comparison
alleviates cognitive burden and elicits respondents’ preferences more
accurately than a choice from a large action space (a demonstration is in Fig.
1). $Z_{1}(u)$ and $Z_{2}(u)$ differ in feature $k$ with probability
$1-|2\pi_{k}-1|$, which increases monotonically from 0 to 1 with the variance
$\pi_{k}(1-\pi_{k})$. It is aligned with the intuition that in order to
maximize the information gain from each question, features with high
uncertainty should have a high chance of being asked. Though the quantity
$p(\mathcal{A}_{i})$ is unknown, it is a constant shared by all the elements
of the stochastic gradient $\tilde{g}$. Therefore, the constant can be
absorbed in the stepsize of the SGD and does not affect the convergence.
Finally, we propose the gradient estimate for GBS with $Z_{1}(u),Z_{2}(u)$
defined in Lemma 1 as

$\displaystyle
g_{\text{GBS}}=(2Y_{i}(Z_{1}(u),Z_{2}(u))-1)(u-\frac{1}{2}),\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
u\sim\prod_{k=1}^{K}\text{Unif}(0,1).$ (3)

The steps of GBS are summarized in Alg. 1.

The data collection of GBS shares the same form of survey with the paired CBC,
which has been widely applied to understanding individual preferences in
marketing, politics, and other computational social science for decades (Luce
and Tukey, , 1964; Toubia et al., , 2003; Hainmueller and Hopkins, , 2015;
Egami and Imai, , 2019; Goplerud et al., , 2022). GBS is thus compatible with
many existing survey systems. The paired choice design used in practice
requires a single choice for each question and rules out the situation of
choosing both or none, which echos the step of zero gradients marginalization
in Eq. 2. In a nutshell, GBS uses the gradient information to automatically
and adaptively design experiments, and in the meanwhile, uses the data from
the experiments to estimate a low-variance stochastic gradient for the optimal
product identification.

input : Number of features $K$, number of questions per respondent $n_{q}$,
stepsize $\eta$

Initialize the logits $\phi$ randomly

while  _not converged_ do

Sample a random individual $i$ from the population.

for  _$j=1,2\cdots n_{q}$_ do

Sample $u\sim\prod_{k=1}^{K}\text{Unif}(0,1)$

Generate a pair of products
$Z_{1}(u)=\mathbf{1}[u>\sigma(-\phi)],Z_{2}(u)=\mathbf{1}[u<\sigma(\phi)]$

Record the respondent’s choice $Y_{i}(Z_{1}(u),Z_{2}(u))$

Compute the gradient estimate
$g_{\text{GBS}}=(2Y_{i}(Z_{1}(u),Z_{2}(u))-1)(u-\frac{1}{2})$

Update $\phi\leftarrow\phi+\eta g_{\text{GBS}}$

end for

end while

Algorithm 1 Gradient-based Survey for Product Design

##  4 Individualized Policy Learning

A single product is often not optimal for all users. For example, a
personalized advertisement email designed based on individual shopping history
may improve the consumers’ open rate. We consider the problem of learning an
individualized policy that assigns a customized product to each user.

Suppose the covariates $X_{i}$ of individual $i$ are observed. The
optimization objective is

$\displaystyle\max_{\theta}V(\theta)=\mathbb{E}_{X_{i}\sim
p(X)}\mathbb{E}_{Z_{i}\sim
p(Z;\pi_{\theta}(X_{i}))}\mathbb{E}[Y_{i}(Z_{i},Z_{0})],$ (4)

where
$p(Z;\pi_{\theta}(X_{i}))=\prod_{k=1}^{K}\text{Bern}(z_{k};\pi_{\theta}(X_{i})_{k})$.
The policy $\pi_{\theta}(X_{i})=\sigma(\phi_{i})$ and the logits
$\phi_{i}=g(X_{i};\theta)\in\mathbb{R}^{K}$ are the output of an amortized
neural network parameterized by $\theta$ with input $X_{i}$.

Applying the results in § 3, an unbiased Monte Carlo estimate of the gradient
w.r.t. the logits $\phi_{i}$ is
$c_{i}(2Y_{i}(Z_{1}(u),Z_{2}(u))-1)(u-\frac{1}{2})$ where
$u\sim\prod_{k=1}^{K}\text{Unif}(0,1)$ and $c_{i}$ is a scalar to be absorbed
in the stepsize. Using the chain rule, the GBS gradient of policy parameter
$\theta$ is

$\displaystyle
g_{\text{GBS}}(\theta)=(2Y_{i}(Z_{1}(u),Z_{2}(u))-1)(u-\frac{1}{2})^{\top}\frac{\partial
g(X_{i};\theta)}{\partial\theta},\quad u\sim\prod_{k=1}^{K}\text{Unif}(0,1)$
(5)

with $Z_{1}(u),Z_{2}(u)$ defined in Lemma 1. When the action space is
combinatorial, existing policy learning methods using offline collected data
often face the extrapolation problem because there is not enough variation in
the observed sets of items (Sachdeva et al., , 2020). We are the first to
enable policy learning with combinatorial actions using conjoint experiments.
The algorithm is summarized in Alg. 2 in Appendix B.

##  5 Related Work

### Conjoint analysis.

Conjoint analysis (CA) uses a survey-based experiment design to measure
multidimensional preferences (Luce and Tukey, , 1964; Toubia et al., , 2003).
We collect the data with a paired comparison survey similar to CA. Different
from CA, GBS does not assume a linear additive utility model and is applicable
to cases with nonlinear interactive utilities. CA estimates partworths as the
choice model parameters (Hainmueller and Hopkins, , 2015) while GBS does not
estimate a choice model and identifies the optimal product using the choice
data directly.

### Adaptive experimental design.

Adaptive conjoint analysis (ACA) progressively refines the attribute levels
presented to respondents for more accurate and efficient data collection
(Toubia et al., , 2004). For example, D-efficiency is designed to maximize
Fisher information (Kuhfeld et al., , 1994), polyhedral methods combine
geometric intuition and analytic center technique to shrink feasible region
(Toubia et al., , 2003, 2007; Sauré and Vielma, , 2019), and adaptive self-
explication integrates attribute importance in design (Netzer and Srinivasan,
, 2011). It might be challenging to generalize heuristics to high-dimensional
features. GBS design is derived from variance reduction of score function
gradient and is aligned with the uncertainty reduction intuitions. Apart from
ACA, an active learning approach is proposed to learn nonparametric choice
models using a directed acyclic graph. Nodes in such a graph are alternative
profiles, which limit the node number to a small scale (Susan et al., , 2022).

### Product design by optimization.

Existing product design methods often adopt optimization heuristics like
Genetic Algorithms (Balakrishnan and Jacob, , 1996; Balakrishnan et al., ,
2004), simulated annealing (Tsafarakis, , 2016), evolutionary algorithm, and
beam search (Paetz et al., , 2021; Hauser, , 2011). Though achieving empirical
improvements, the properties and generalization abilities of the heuristics is
largely unclear. Discrete optimization methods such as Lagrangian relaxation
with branch-and-bound have also been applied (Camm et al., , 2006; Belloni et
al., , 2008). More broadly, deep learning methods such as variational
autoencoders are used for the design of product aesthetics (Burnap et al., ,
2023).

### Policy learning.

The customized product design is related to policy learning. The offline
policy learning often estimates the outcome function (Wang et al., , 2016) or
directly optimizes the value function by propensity weighting (Athey and
Wager, , 2017). However, it is hardly possible for the observational data to
contain the outcomes of all the actions from a combinatorial space. The
learned policy might be suboptimal due to the overlapping issue. Online policy
learning uses bandits and reinforcement learning to maximize the cumulative
reward. Qin et al.,  (2014) explores combinatorial action spaces, which
require users to select actions from the complete action space and need
external covariates for each feature dimension. GSB overcomes these challenges
using the partial profile design from the conjoint analysis.

##  6 Empirical Study

We compare GBS with baseline models using simulated data similar to Toubia et
al.,  (2007). The product $Z$ is represented as $K$ binary features. First, we
study the problem of identifying a single optimal product. For each
independent trial, the population-level marginal preferences for the features
are generated by $\mu\sim\mathcal{N}(a,I)$, $a=(1,\cdots,1)$. The individual
preferences (partworths) are generated by $W_{i}\sim\mathcal{N}(\mu,I)$ for
each individual $i$.

### Data generation.

We consider three types of representative utilities $V_{i}(Z)$. The first type
is linear utility $V_{i}(Z)=W_{i}^{\top}Z$. The second type includes pairwise
interactions
$V_{i}(Z)=W_{i}^{\top}Z+\sum_{k,k^{\prime}}\tilde{W}_{i}^{kk^{\prime}}Z_{k}Z_{k^{\prime}}$
where $\tilde{W}_{i}^{kk^{\prime}}\sim\text{Unif}(-2a,0)$. For $K=10$, all the
pairwise interactions are included, and for $K=100$, a subset of 100 pairwise
interaction terms is used to compute $V_{i}(Z)$. The third type set
$V_{i}(Z)=f_{0}(Z)$ where $f_{0}(\cdot)$ is a pre-trained neural network. This
type of utility includes higher-order interactions of the product features.
The data for each method are the query product features and the choices
collected from paired choice questions, i.e.,
$\\{(Y_{i}(Z_{1}^{ij},Z_{0}^{ij}),Z_{1}^{ij},Z_{0}^{ij})\\}_{i=1:N}^{j=1:n_{q}}$.
For non-adaptive methods, the item pair in a question is generated randomly.
Each respondent makes $n_{q}=10$ times the choices.

### Baselines.

A logistic model (Logistic) assumes $\hat{V}_{i}(Z)=W^{\top}Z$. For a pair of
products $(Z_{1},Z_{2})$, the likelihood of choosing $Z_{1}$ is
$1/(1+\exp(W^{\top}(Z_{2}-Z_{1})))$. The parameter $W$ is estimated by maximum
likelihood and the optimal product is $\mathbf{1}[\hat{W}_{MLE}>0]$. A mixed
logit model is a hierarchical Bayes (HB) model widely used in conjoint
analysis (Allenby et al., , 2005). It assumes $m\sim\mathcal{N}(0,I)$,
$w_{i}\,|\,m\sim\mathcal{N}(\mu,I)$, and
$p(Y(Z_{1},Z_{2})=1\,|\,w_{i},Z_{1},Z_{2})=1/(1+\exp(w_{i}^{\top}(Z_{2}-Z_{1})))$.
We estimate $m$ as the maximum a posteriori estimation using the PyMC package,
and the estimated optimal product is $\mathbf{1}[\hat{m}_{MAP}>0]$. Another
baseline takes the representative utility $\hat{V}_{i}(Z)=f_{\gamma}(Z)$.
$f_{\gamma}(\cdot)$ is a feedforward neural network (NN) with two hidden
layers and parameters $\gamma$. $\gamma$ is estimated as an MLE. The estimated
optimal product is $\operatorname*{arg\,max}_{Z}f_{\gamma}(Z)$, which requires
enumerating all possible products.

### Evaluation metrics.

The first metric for evaluating a chosen product is the average utility on a
hold-out test set. It compares the relative performance of different methods.
When computationally feasible, we also rank all the possible products
according to their average utility on the population from high to low. The
ranking allows the evaluation of absolute performance compared to the global
optimum.

(a) Type 1 (b) Type 2 (c) Type 3 Figure 2: The utility (Top row, higher is
better) and the ranking (Bottom row, lower is better) for the selected
products with $K=10$ features. GBS has similar performance with correctly
specified choice utility models (Type 1), and identifies better product when
the true utility is unknown (Type 2,3).

§ 6 shows the test utility and ranking of the estimated optimal product across
a different number of respondents for a product with $K=10$ features. When the
utility function is correctly specified (Type 1), Logistic and HB reach the
highest utility with a small number of respondents. However, when the utility
model is misspecified (Type 2, 3), the performance of Logistic and HB with
linearity assumption significantly deteriorates. For Type 2 utility, the
ranking drops with an increasing number of respondents when the driving factor
for the estimates shifts from variance to bias. NN assumes a flexible
nonparametric utility function, but when the sample size is small, its
performance is dominated by the variances in the data. Fig. 5 in the Appendix
B shows the results of NN with additional respondents. For the utility
function without interactions, it needs 500 respondents to select the optimal
product, and for relatively complex Type 2 and 3 utility, the number of
respondents needed is 4000 and 2000, respectively. The weak data efficiency
increases the experimental cost and may be infeasible. Moreover, finding the
optimal product with a trained NN utility function needs to explore all
possible products as inputs, which becomes challenging when the feature
dimension is high. In comparison, GBS identifies the optimal products across
all utility types with less than 100 respondents. For the linear utility, when
the respondent number is small, the correctly specified models outperform GBS,
but the performance gap diminishes quickly when the respondent number
increases to around 70. For nonlinear utilities, GBS outperforms the baseline
models and achieves the global optimum.

(a) Type 1

(b) Type 2

(c) Type 3

Figure 3: The utility (higher is better) with $K=100$ features.

Figure 4: Personalized Product.

Fig. 4 shows the results for the single product design with $K=100$ features.
Finding the optimal product with a trained NN needs to compare the predicted
utility of $2^{100}$ items, which is infeasible to compute, so we drop it from
the baselines. Similarly, we drop the ranking metric that needs to evaluate
all product combinations. The performance pattern is similar to $K=10$. GBS
can identify the optimal product from around $10^{30}$ choices in less than
half a minute. It is flexible with the underlying utility function and is
efficient with data size, which makes it practically applicable.

§§ 6 and 4 illustrates a trade-off between specification and estimation (Susan
et al., , 2022). A complex utility model might capture a broader spectrum of
choice behaviors and hence has a small specification error; it could also
introduce big estimation errors when limited to a small data set. Instead of
estimating a full choice model, GBS directly uses the choice data for product
design, bypassing this trade-off. GBS is flexible with the underlying utility
function and is efficient with data size, which makes it practically
applicable.

Next, we study the personalized product design with choice experiments. The
data is generated the same as before, except the individual preferences follow
a mixture distribution $W_{i}\sim
0.5\mathcal{N}(\mu_{1},\Sigma)+0.5\mathcal{N}(\mu_{2},\Sigma)$ where $\mu_{1}$
is positive on the first half of elements and negative on the others, and
$\mu_{2}$ is opposite. The mixture distribution reflects the preference
heterogeneity in the population. We assume the observed covariates
$X_{i}=\exp(W_{i})$ of individual $i$ and take the utility as the nonlinear
Type 2. We modify the neural network utility model as
$\hat{V}(Z,X_{i})=Z^{\top}f_{\gamma}(X_{i})$ (denoted as NN-ind). The utility
form gives an estimated optimal product for individual $i$ as
$\mathbf{1}[f_{\gamma}(X_{i})>0]$ without the need to evaluate all possible
products’ utility for each individual.

Fig. 4 shows the selected products’ utility averaged on a hold-out test set.
Since Logistic estimates a single optimal product, the estimate is a product
with zero utility due to the symmetry of the population. NN-ind adapts to
individual heterogeneity and has utility higher than Logistic. However, it is
subject to model misspecification due to the last linear layer and may have
low data efficiency as in the case of the single optimal product. GBS reaches
the highest utility, and the performance improves with more respondents
joining the experiments.

##  7 Discussions

This paper bridges the domains of gradient-based machine learning and discrete
choice experiments. The proposed GBS is flexible with the underlying form of
choice utility, is data-efficient by an adaptive design, is scalable to high-
dimensional features, and is applicable to uniform or personalized product
designs.

However, there is no free lunch. GBS does not estimate a choice model. It does
not provide a full rank of all the possible products. To explain people’s
preferences, GBS may need explainable AI techniques such as saliency map
(Adebayo et al., , 2018) rather than estimating the preferences as a part of
model parameters. Nevertheless, GBS provides a flexible optimization
framework. Except for maximizing the market share, the objective may
incorporate the costs and prices of a product to maximize the profit. If the
product design is under a budget constraint, the constrained optimization
might be considered using methods such as projected gradient.

From a manager’s view, a company often needs to design a product line
consisting of several products (Balakrishnan et al., , 2004; Belloni et al., ,
2008). One way to apply GBS for this task is by a separate approach, where the
population is clustered into segments, and a single best product is determined
for each segment (Paetz et al., , 2021). It is also feasible to model a
product line as several binary vectors and apply GBS to design a product line
jointly.

GBS builds on the inference of discrete latent variables in machine learning
(Gu et al., , 2015; Jang et al., , 2017; Tucker et al., , 2017; Kool et al., ,
2019; Dong et al., , 2020; Kunes et al., , 2023). If a gradient estimator
contrasts several values of the objective function, it can potentially be used
for the adaptive question design. Combining recent discrete optimization
techniques with discrete choice experiments is an interesting future
direction.

## References

  * Adebayo et al., (2018) Adebayo, J., Gilmer, J., Muelly, M., Goodfellow, I., Hardt, M., and Kim, B. (2018).  Sanity checks for saliency maps.  Advances in neural information processing systems, 31.
  * Allenby et al., (2005) Allenby, G. M., Rossi, P. E., and McCulloch, R. E. (2005).  Hierarchical bayes models: A practitioners guide.  SSRN Electronic Journal.
  * Athey and Wager, (2017) Athey, S. and Wager, S. (2017).  Policy learning with observational data.
  * Balakrishnan et al., (2004) Balakrishnan, P., Gupta, R., and Jacob, V. S. (2004).  Development of hybrid genetic algorithms for product line designs.  IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics), 34(1):468–483.
  * Balakrishnan and Jacob, (1996) Balakrishnan, P. and Jacob, V. S. (1996).  Genetic algorithms for product design.  Management Science, 42(8):1105–1117.
  * Belloni et al., (2008) Belloni, A., Freund, R., Selove, M., and Simester, D. (2008).  Optimizing product line designs: Efficient methods and comparisons.  Management Science, 54(9):1544–1552.
  * Bertsimas and Mišić, (2017) Bertsimas, D. and Mišić, V. V. (2017).  Robust product line design.  Operations Research, 65(1):19–37.
  * Burnap et al., (2023) Burnap, A., Hauser, J. R., and Timoshenko, A. (2023).  Product aesthetic design: A machine learning augmentation.  Marketing Science.
  * Camm et al., (2006) Camm, J. D., Cochran, J. J., Curry, D. J., and Kannan, S. (2006).  Conjoint optimization: An exact branch-and-bound algorithm for the share-of-choice problem.  Management Science, 52(3):435–447.
  * Dong et al., (2020) Dong, Z., Mnih, A., and Tucker, G. (2020).  DisARM: an antithetic gradient estimator for binary latent variables.
  * Egami and Imai, (2019) Egami, N. and Imai, K. (2019).  Causal interaction in factorial experiments: Application to conjoint analysis.  Journal of the American Statistical Association, 114(526):529–540.
  * Ellickson et al., (2022) Ellickson, P. B., Kar, W., and Reeder, J. C. (2022).  Estimating marketing component effects: Double machine learning from targeted digital promotions.  Marketing Science.
  * Goplerud et al., (2022) Goplerud, M., Imai, K., and Pashley, N. E. (2022).  Estimating heterogeneous causal effects of High-Dimensional treatments: Application to conjoint analysis.  arXiv preprint arXiv:2201.01357.
  * Green et al., (1991) Green, P. E., Krieger, A. M., and Agarwal, M. K. (1991).  Adaptive conjoint analysis: Some caveats and suggestions.  Journal of Marketing Research, 28(2):215–222.
  * Green and Rao, (1971) Green, P. E. and Rao, V. R. (1971).  Conjoint measurement-for quantifying judgmental data.  Journal of Marketing research, 8(3):355–363.
  * Green and Srinivasan, (1978) Green, P. E. and Srinivasan, V. (1978).  Conjoint analysis in consumer research: issues and outlook.  Journal of consumer research, 5(2):103–123.
  * Gu et al., (2015) Gu, S., Levine, S., Sutskever, I., and Mnih, A. (2015).  Muprop: Unbiased backpropagation for stochastic neural networks.  arXiv preprint arXiv:1511.05176.
  * Hainmueller and Hopkins, (2015) Hainmueller, J. and Hopkins, D. J. (2015).  The hidden american immigration consensus: A conjoint analysis of attitudes toward immigrants.  American journal of political science, 59(3):529–548.
  * Hauser, (2011) Hauser, J. R. (2011).  New developments in product-line optimization.
  * Hein et al., (2020) Hein, M., Kurz, P., and Steiner, W. J. (2020).  Analyzing the capabilities of the hb logit model for choice-based conjoint analysis: a simulation study.  Journal of Business Economics, 90(1):1–36.
  * Huber, (2005) Huber, J. (2005).  Conjoint analysis: how we got here and where we are (an update).  In Sawtooth software conference, volume 98382, page 360. Citeseer.
  * Jang et al., (2017) Jang, E., Gu, S., and Poole, B. (2017).  Categorical reparameterization with Gumbel-softmax.  In 5th International Conference on Learning Representations, ICLR 2017, Toulon, France, April 24-26, 2017, Conference Track Proceedings. OpenReview.net.
  * Kohli and Krishnamurti, (1989) Kohli, R. and Krishnamurti, R. (1989).  Optimal product design using conjoint analysis: Computational complexity and algorithms.  European Journal of Operational Research, 40(2):186–195.
  * Kool et al., (2019) Kool, W., van Hoof, H., and Welling, M. (2019).  Buy 4 reinforce samples, get a baseline for free!
  * Kuhfeld et al., (1994) Kuhfeld, W. F., Tobias, R. D., and Garratt, M. (1994).  Efficient experimental design with marketing research applications.  Journal of Marketing Research, 31(4):545–557.
  * Kunes et al., (2023) Kunes, R. Z., Yin, M., Land, M., Haviv, D., Pe’er, D., and Tavaré, S. (2023).  Gradient estimation for binary latent variables via gradient variance clipping.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 37, pages 8405–8412.
  * Louviere and Woodworth, (1983) Louviere, J. J. and Woodworth, G. (1983).  Design and analysis of simulated consumer choice or allocation experiments: an approach based on aggregate data.  Journal of marketing research, 20(4):350–367.
  * Luce and Tukey, (1964) Luce, R. D. and Tukey, J. W. (1964).  Simultaneous conjoint measurement: A new type of fundamental measurement.  Journal of mathematical psychology, 1(1):1–27.
  * McFadden, (1974) McFadden, D. (1974).  The measurement of urban travel demand.  Journal of public economics, 3(4):303–328.
  * Netzer and Srinivasan, (2011) Netzer, O. and Srinivasan, V. (2011).  Adaptive self-explication of multiattribute preferences.  Journal of Marketing Research, 48(1):140–156.
  * Paetz et al., (2021) Paetz, F., Steiner, W., and Hruschka, H. (2021).  Optimal product line designing: A comparison of different discrete choice-based approaches.  In EMAC 2021 Annual Conference.
  * Qin et al., (2014) Qin, L., Chen, S., and Zhu, X. (2014).  Contextual combinatorial bandit and its application on diversified online recommendation.  In Proceedings of the 2014 SIAM International Conference on Data Mining, pages 461–469. SIAM.
  * Sachdeva et al., (2020) Sachdeva, N., Su, Y., and Joachims, T. (2020).  Off-policy bandits with deficient support.  In Proceedings of the 26th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining, pages 965–975.
  * Sauré and Vielma, (2019) Sauré, D. and Vielma, J. P. (2019).  Ellipsoidal methods for adaptive choice-based conjoint analysis.  Operations Research, 67(2):315–338.
  * Susan et al., (2022) Susan, F., Golrezaei, N., Emamjomeh-Zadeh, E., and Kempe, D. (2022).  Active learning for non-parametric choice models.  arXiv preprint arXiv:2208.03346.
  * Toubia et al., (2007) Toubia, O., Hauser, J., and Garcia, R. (2007).  Probabilistic polyhedral methods for adaptive choice-based conjoint analysis: Theory and application.  Marketing Science, 26(5):596–610.
  * Toubia et al., (2004) Toubia, O., Hauser, J. R., and Simester, D. I. (2004).  Polyhedral methods for adaptive choice-based conjoint analysis.  Journal of Marketing Research, 41(1):116–131.
  * Toubia et al., (2003) Toubia, O., Simester, D. I., Hauser, J. R., and Dahan, E. (2003).  Fast polyhedral adaptive conjoint estimation.  Marketing Science, 22(3):273–303.
  * Train, (2009) Train, K. E. (2009).  Discrete choice methods with simulation.  Cambridge university press.
  * Tsafarakis, (2016) Tsafarakis, S. (2016).  Redesigning product lines in a period of economic crisis: a hybrid simulated annealing algorithm with crossover.  Annals of Operations Research, 247:617–633.
  * Tucker et al., (2017) Tucker, G., Mnih, A., Maddison, C. J., Lawson, J., and Sohl-Dickstein, J. (2017).  Rebar: Low-variance, unbiased gradient estimates for discrete latent variable models.  In Advances in Neural Information Processing Systems, pages 2627–2636.
  * Wang et al., (2016) Wang, Y., Wu, P., Liu, Y., Weng, C., and Zeng, D. (2016).  Learning optimal individualized treatment rules from electronic health record data.  IEEE International Conference on Healthcare Informatics. IEEE International Conference on Healthcare Informatics, 2016:65–71.
  * Yin et al., (2020) Yin, M., Ho, N., Yan, B., Qian, X., and Zhou, M. (2020).  Probabilistic best subset selection via gradient-based optimization.  arXiv preprint arXiv:2006.06448.
  * Yin and Zhou, (2019) Yin, M. and Zhou, M. (2019).  Arm: Augment-reinforce-merge gradient for stochastic binary networks.  In International Conference on Learning Representations.

##  Appendix A Proof

This section provides the proof of Lemmas 1, 2 and 3.

###### Proof of Lemma 1.

The optimization objective is

$\displaystyle V(\phi)=\mathbb{E}_{Z\sim
p(Z;\sigma(\phi))}\mathbb{E}[Y_{i}(Z,Z_{0})\,|\,Z],$

where $p(Z;\sigma(\phi))=\prod_{k=1}^{K}\text{Bern}(z_{k};\sigma(\phi_{k}))$,
$Z=(z_{1},\cdots,z_{K})$. The k-th element for the gradient of $V(\phi)$ is

$\displaystyle\nabla_{\phi_{k}}V(\phi)$
$\displaystyle=\mathbb{E}_{\boldsymbol{z}_{\backslash k}\sim\prod_{\nu\neq
k}{\text{Bern}}(z_{\nu};\sigma(\phi_{\nu}))}\\{\nabla_{\phi_{k}}\mathbb{E}_{z_{k}\sim{\text{Bern}}(\sigma(\phi_{k}))}\mathbb{E}[Y_{i}(Z,Z_{0})\,|\,Z]\\}$
(6)

Denote $f(z_{k})=\mathbb{E}[Y_{i}(Z,Z_{0})\,|\,Z]$, we have

$\displaystyle\nabla_{\phi_{k}}\mathbb{E}_{z_{k}\sim{\text{Bern}}(\sigma(\phi_{k}))}[f(z_{k})]$
$\displaystyle=$ $\displaystyle\sigma(\phi_{k})\sigma(-\phi_{k})[f(1)-f(0)]$
$\displaystyle=$
$\displaystyle\mathbb{E}_{u\sim\text{Unif}(0,1)}[f(\mathbf{1}{[u<\sigma(\phi)]})(1-2u)]$
$\displaystyle=$
$\displaystyle\mathbb{E}_{u\sim\text{Unif}(0,1)}[f(\mathbf{1}{[u<\sigma(\phi)]})(1/2-u)]+\mathbb{E}_{\tilde{u}\sim\text{Unif}(0,1)}[f(\mathbf{1}{[\tilde{u}<\sigma(\phi)]})(1/2-\tilde{u})]$
$\displaystyle=$
$\displaystyle\mathbb{E}_{u\sim\text{Unif}(0,1)}\left[\big{(}f(\mathbf{1}{[u>\sigma(-\phi)]})-f(\mathbf{1}{[u<\sigma(\phi)]})\big{)}(u-1/2)\right].$
(7)

The first two equations can be straightfowardly evaluated since it is an
expectation with a scalar variable. The third equation is applying antithetic
sampling with $\tilde{u}=1-u$. The last equation is summing up the two
expectations. Plugging Eq. 7 into Eq. 6 gives

$\displaystyle\nabla_{\phi_{k}}V(\phi)$ $\displaystyle=$
$\displaystyle\mathbb{E}_{\boldsymbol{z}_{\backslash k}\sim\prod_{\nu\neq
k}{\text{Bern}}(z_{\nu};\sigma(\phi_{\nu}))}\big{\\{}\mathbb{E}_{u_{k}\sim{\text{Unif}}(0,1)}\Big{[}\big{(}\mathbb{E}[Y_{i}(Z,Z_{0})\,|\,Z=(\boldsymbol{z}_{\backslash
k},z_{k}=\mathbf{1}{[u_{k}>\sigma(-\phi_{k})]})]-$
$\displaystyle\mathbb{E}[Y_{i}(Z,Z_{0})\,|\,Z=(\boldsymbol{z}_{\backslash
k},z_{k}=\mathbf{1}{[u_{k}<\sigma(\phi_{k})]})]\big{)}(u_{k}-\frac{1}{2})\Big{]}$
$\displaystyle=$
$\displaystyle\mathbb{E}_{\boldsymbol{u}\sim\prod_{k=1}^{K}{\text{Unif}}(0,1)}\Big{[}\big{(}\mathbb{E}[Y_{i}(Z,Z_{0})\,|\,Z=\mathbf{1}{[\boldsymbol{u}>\sigma(-\boldsymbol{\phi})]})]-\mathbb{E}[Y_{i}(Z,Z_{0})\,|\,Z=\mathbf{1}{[\boldsymbol{u}<\sigma(\boldsymbol{\phi})]})]\big{)}(u_{k}-\frac{1}{2})\Big{]}$

Therefore,

$\displaystyle\nabla_{\phi}V(\phi)=\mathbb{E}_{u\sim\prod_{k=1}^{K}\text{Unif}(0,1)}\left[\mathbb{E}[(Y_{i}(Z_{1}(u),Z_{0})-Y_{i}(Z_{2}(u),Z_{0}))(u-\frac{1}{2})\,|\,u]\right],$
(8)

where
$Z_{1}(u)=\mathbf{1}[u>\sigma(-\phi)],Z_{2}(u)=\mathbf{1}[u<\sigma(\phi)].$ ∎

###### Proof of Lemma 2.

Denote $S=Y(Z_{1},Z_{0})-Y(Z_{2},Z_{0})$, $T=2Y(Z_{1},Z_{2})-1$. We have

$\displaystyle p(S=1)=$ $\displaystyle p(Y(Z_{1},Z_{0})=1,Y(Z_{2},Z_{0})=0)$
$\displaystyle=$
$\displaystyle\frac{e^{U(Z_{1})}}{e^{U(Z_{1})}+e^{U(Z_{0})}}\cdot\frac{e^{U(Z_{0})}}{e^{U(Z_{2})}+e^{U(Z_{0})}}$

and

$\displaystyle p(S=1\,|\,Y(Z_{1},Z_{0})\neq Y(Z_{2},Z_{0}))$ $\displaystyle=$
$\displaystyle\frac{p(Y(Z_{1},Z_{0})=1,Y(Z_{2},Z_{0})=0)}{p(Y(Z_{1},Z_{0})=1,Y(Z_{2},Z_{0})=0)+p(Y(Z_{1},Z_{0})=0,Y(Z_{2},Z_{0})=1)}$
$\displaystyle=$
$\displaystyle\frac{e^{U(Z_{1})+U(Z_{0})}}{e^{U(Z_{1})+U(Z_{0})}+e^{U(Z_{2})+U(Z_{0})}}$
$\displaystyle=$ $\displaystyle\frac{e^{U(Z_{1})}}{e^{U(Z_{1})}+e^{U(Z_{2})}}$
$\displaystyle=$ $\displaystyle p(Y(Z_{1},Z_{2})=1)$ $\displaystyle=$
$\displaystyle p(T=1).$

Similarly, we have $p(S=-1\,|\,Y(Z_{1},Z_{0})\neq Y(Z_{2},Z_{0}))=p(T=-1)$.
Therefore, $S\,|\,Y(Z_{1},Z_{0})\neq Y(Z_{2},Z_{0})\stackrel{{\scriptstyle
d}}{{=}}T$. ∎

input : Individual covariates $\\{X_{i}\\}$, policy function
$g(\cdot;\theta)$, number of features $K$, number of questions per respondent
$n_{q}$, stepsize $\eta$

Initialize the policy parameters $\theta$ randomly.

while  _not converged_ do

Sample a random individual $i$ from the population.

for  _$j=1,2\cdots n_{q}$_ do

Sample $u\sim\prod_{k=1}^{K}\text{Unif}(0,1)$

Generate a pair of products
$Z_{1}(u)=\mathbf{1}[u>\sigma(-\phi)],Z_{2}(u)=\mathbf{1}[u<\sigma(\phi)]$

Record the respondent’s choice $Y_{i}(Z_{1}(u),Z_{2}(u))$

Compute the gradient estimate $g_{\text{GBS}}(\theta)$ by Eq. 11

Update $\theta\leftarrow\theta+\eta g_{\text{GBS}}(\theta)$

end for

end while

Algorithm 2 Gradient-based Survey for Policy Learning

###### Proof of Lemma 3.

Denote the event $\mathcal{A}_{i}=\\{Y_{i}(Z_{1},Z_{0})\neq
Y_{i}(Z_{2},Z_{0})\\}$, we have

$\displaystyle\mathbb{E}_{u\sim\prod_{k=1}^{K}\text{Unif}(0,1)}\left[\mathbb{E}\Big{[}(Y_{i}(Z_{1}(u),Z_{0})-Y_{i}(Z_{2}(u),Z_{0}))(u-\frac{1}{2})\,|\,u\Big{]}\right]$
(9) $\displaystyle=$
$\displaystyle\mathbb{E}_{u\sim\prod_{k=1}^{K}\text{Unif}(0,1)}\left[\mathbb{E}\Big{[}(Y_{i}(Z_{1}(u),Z_{0})-Y_{i}(Z_{2}(u),Z_{0}))(u-\frac{1}{2})\,|\,u,\mathcal{A}_{i}\Big{]}p(\mathcal{A}_{i})\right]$
(10) $\displaystyle=$
$\displaystyle\mathbb{E}_{u\sim\prod_{k=1}^{K}\text{Unif}(0,1)}\left[\mathbb{E}\Big{[}(2Y_{i}(Z_{1}(u),Z_{2}(u))-1)(u-\frac{1}{2})p(\mathcal{A}_{i})\,|\,u\Big{]}\right].$
(11)

The first equality is by the law of total expectation; the second is by Lemma
2 and the law of unconscious statistician (LOTUS). The Monte Carlo estimation
of the gradient in Eq. 11 is
$\tilde{g}=(2Y_{i}(Z_{1}(u),Z_{2}(u))-1)(u-\frac{1}{2})p(\mathcal{A}_{i})$,
$u\sim\prod_{k=1}^{K}\text{Unif}(0,1)$. ∎

##  Appendix B Additional Results

This section contains the customized product design algorithm and additional
results for the empirical study in § 6. Fig. 5 shows the utility and ranking
for the NN baseline with a large number of respondents.

(a) Type 1

(b) Type 2

(c) Type 3

Figure 5: The utility (Top row) and the ranking (Bottom row) for the NN
baseline with a large number of respondents.