# Consistent Training via Energy-Based GFlowNets for Modeling Discrete Joint
Distributions

Chanakya Ekbote1, Moksh Jain1, Payel Das2, Yoshua Bengio1,3
1Mila, Université de Montréal 2IBM Thomas J Watson Research Center 3 CIFAR AI
Chair Correspondence<EMAIL_ADDRESS>
###### Abstract

Generative Flow Networks (GFlowNets) have demonstrated significant performance
improvements for generating diverse discrete objects $x$ given a reward
function $R(x)$, indicating the utility of the object and trained
independently from the GFlowNet by supervised learning to predict a desirable
property $y$ given $x$. We hypothesize that this can lead to incompatibility
between the inductive optimization biases in training $R$ and in training the
GFlowNet, potentially leading to worse samples and slow adaptation to changes
in the distribution. In this work, we build upon recent work on jointly
learning energy-based models with GFlowNets and extend it to learn the joint
over multiple variables, which we call Joint Energy-Based GFlowNets (JEBGFNs),
such as peptide sequences and their antimicrobial activity. Joint learning of
the energy-based model, used as a reward for the GFlowNet, can resolve the
issues of incompatibility since both the reward function $R$ and the GFlowNet
sampler are trained jointly. We find that this joint training or joint energy-
based formulation leads to significant improvements in generating anti-
microbial peptides. As the training sequences arose out of evolutionary or
artificial selection for high antibiotic activity, there is presumably some
structure in the distribution of sequences that reveals information about the
antibiotic activity. This results in an advantage to modeling their joint
generatively vs. pure discriminative modeling. We also evaluate JEBGFN in an
active learning setting for discovering anti-microbial peptides.

##  1 Introduction

Biological sequence design is a problem of paramount importance for various
use cases such as design of new materials, discovering new drugs, and so on.
However, designing these sequences is a time consuming process requiring a lot
of domain expertise and many expensive wet lab experiments. Depending on the
use case, this can involve searching over an exponentially large combinatorial
space (upto orders of $20^{60}$).

Over the past couple of years, a spate of works [1, 7, 19, 8, 27, 26] have
leveraged recent progress in deep generative methods to aid the search for
viable biological sequences. These methods, which fall primarily in the broad
areas of reinforcement learning (RL) and Bayesian optimization (BO), typically
involve finding sequences that _maximize_ certain properties. Consequently,
they suffer from one key drawback, the generated sequences have low diversity,
i.e. the generated sequences are “similar". Considering the diverse nature of
viable biological targets, and the drug or material discovery process,
diversity of generated sequences is a key consideration [20, 12].

Generative Flow Networks [GFlowNets; 2, 3] are a family of probabilistic
models that alleviate this issue by sampling sequences with probability
proportional to the reward, as opposed to reward maximization, resulting in
generation of diverse candidates. GFlowNet-AL [12] leveraged GFlowNets in an
active learning setting to generate _diverse and novel_ biological sequences.
In each round of active learning GFlowNet-AL learns a predictor for estimating
the reward $R$ from the oracle using the data collected so far. This predictor
is then used as the reward for training a GFlowNet. Since the reward function
and the GFlowNet are learned separately, we hypothesize that this can lead to
incompatibility between the reward $R$ and the GFlowNet, potentially leading
to poor quality of samples being generated and slow adaptation to changes in
the distribution. Moreover, since the generation of sequences is only based on
a single reward, this method cannot be used to generate sequences with a
specific set of desired properties. For instance, a practical problem of
interest is to generate peptide sequences that target a specific class of
bacteria, rather than having general toxicity.

Energy-based models (EBMs) [15, 10, 24] have recently been successfully
applied on a wide variety of probabilistic modelling and inference tasks. In
particular, [16, 28] demonstrate that EBMs adapt quickly to changes in the
distribution in continual setting. [29] proposed a framework for jointly
training an energy function with a GFlowNet sampler. The GFlowNet samples data
proportional to the reward (exp(-energy)) associated with the learned energy
function, which models the unnormalised data distribution. The energy function
is learnt via contrastive divergence, and the required negative samples are
generated by the GFlowNet. Since the reward $R$ (i.e. the energy) and the
GFlowNet are both trained via an alternating optimization paradigm, we
hypothesize that it reduces the presumed incompatibility faced by training
them separately.

We extend the framework presented in [29] to model the joint distribution over
multiple variables. In particular, we study the problem of learning the joint
$p(x,y)$ where $x$ is the data, and $y$ is the corresponding label. Our work
offers a GFlowNet-based alternative to prior work on modelling the joint [23,
6, 5, 11, 21, 14] for compositional discrete spaces. To summarize, our key
contributions are as follows:

  * •

We introduce Joint Energy-based GFlowNets (JEB-GFN), an extension of Energy-
Based GFlowNets for modelling joint distributions over multiple variables.

  * •

We empirically demonstrate that JEB-GFN can learn the joint distribution over
two variables, data $x$ and its corresponding label $y$, and leverage it for
conditional generation.

  * •

Finally, we evaluate JEB-GFN for generating peptides with anti-microbial
properties in a single round as well as active learning setting.

##  2 Background

Generative Flow Networks [2, 3, GFlowNets;] view the problem of sampling
discrete objects $x\in\mathcal{X}$ as a problems of learning a control policy.
GFlowNets learn a stochastic policy $\pi$ to generate $x$ such that the
marginal likelihood of sampling $x\sim\pi(x)$ is proportional to a non-
negative reward $R(x)$. In this section we review some basic notation and
terminology and refer the reader to [3] for a more thorough discussion.

Generative Flow Networks are a probabilistic framework for learning stochastic
policies for discrete data. They model the discrete generational process of
discrete objects $x\in\mathcal{X}$, although bengio2021gflownet [3] propose a
methodology for also generating continuous-valued objects. A compositional
object $x\in\mathcal{X}$, is modeled as a _trajectory_ comprising of a
sequence of discrete actions from a set of actions $\mathcal{A}$. At each step
of this trajectory, we obtain a partially constructed object, which we call
state $s\in\mathcal{S}$, which is modified by the next action. The space of
possible trajectories is described by a directed acyclic graph (DAG),
$\mathcal{G}=(\mathcal{S},\mathcal{A})$, whose vertices $\mathcal{S}$ are the
set of all possible states (partially constructed objects along with a special
empty state $s_{0}$) and the edges at each state $\mathcal{A}$ are actions
that modify one state to another. Trajectories starting at $s_{0}$ and
terminating at $x\in\mathcal{X}$ are called _complete trajectories_ , denoted
by $\tau\in\mathcal{T}$ where $\mathcal{T}$ is the set of all possible
complete trajectories. Note that while this graph does not need to be
instantiated and can be arbitrary in general. We consider the DAG case, so the
actions can only be constructive.

The _trajectory flow_ is the unnormalized probability density over the set of
all complete trajectories $\tau\in\mathcal{T}$; defined as a non negative
function $F(\tau)$. For some constant $Z$, we can compute the probability
associated with the trajectory as $P_{F}(\tau)=F(\tau)/Z$.
$P_{F}(s_{t+1}|s_{t})$ is called the forward policy, and sampling iteratively
from this forward policy allows one to generate complete trajectories, with a
probability proportional to the unnormalized density $F$. We denote the
marginal probability of sampling complete trajectories following $P_{F}$
terminating in $x$ as $\pi(x)$.

The underlying learning problem solved by GFlowNets is to estimate a $P_{F}$
such that the marginal likelihood of generating object $x$ is proportional to
a given non-negative reward $R(x)$. [2] introduced a flow-matching objective
to learn this policy. In this work, we focus on the Trajectory Balance view
introduced by [17]. We estimate a forward policy $P_{F}$, a backward policy
$P_{B}$ and a partition function $Z$. The backward policy $P_{B}$ is a
distribution over the set of parents given an input state $s$, which
prescribes how to sample trajectories starting at any terminal state $x$ and
terminating in $s_{0}$, i.e. ways in which $x$ _could have been_ constructed.
The partition function $Z$ denotes the total reward
$\sum_{\tau\in\mathcal{T}}F(\tau)=\sum_{x}R(x)$ and is also estimated during
training (by $Z_{\theta}$ below). The parameters $\theta$ parameterizing
$P_{F}(-|s;\theta),P_{B}(-|s;\theta),Z_{\theta}$ are estimated by minimizing
the following per-trajectory loss:

$\mathcal{L}(\tau;\theta)=\left(\log\frac{Z_{\theta}\prod_{t=0}^{t=n-1}P_{F}(s_{t+1}|s_{t};\theta)}{R(s_{n})\prod_{t=0}^{t=n-1}P_{B}(s_{t}|s_{t+1};\theta)}\right)^{2}.$
(1)

[12] introduced GFlowNet-AL, an active learning algorithm leveraging GFlowNets
for generating _diverse_ candidates that maximize some property of the
candidate that is expensive to evaluate (oracle). It consists of several
rounds, each consisting of generating a batch of candidates that are evaluated
with the oracle. To generate each batch of candidates, it first learns a
supervised predictor to approximate the oracle and its epistemic uncertainty
using a dataset of past evaluations of the oracle. This predictor (and the
epistemic uncertainty associated with it) is then converted via an acquisition
function into the reward function for training a GFlowNet. As the GFlowNet and
reward are trained separately, we have observed that adaptation is slow over
rounds (especially if retraining the reward function from scratch) and can
result in inconsistencies between the GFlowNet and the reward. This motivates
the approach proposed here and described below.

##  3 Energy Based Modelling with GFlowNets

As introduced in section 2, the training regime of GFlowNets, requires an
explicit reward function $R(x)$. Leveraging GFlowNets for generative modelling
would require setting the reward function to the unnormalized probability of a
particular sample. In practice, one does not have access to the closed-form
solution of this unnormalised probability, but only to data sampled from a
target distribution. To combat this issue, [29] introduced the framework of
Energy Based GFlowNets (EB-GFNs), where the reward function for the GFlowNet
is set using the output from an energy based model (EBM)
$\mathcal{E}_{\phi}(x)$ learned simultaneously. The final products obtained
include not only be the trained GFlowNet but also the energy based model
$\mathcal{E}_{\phi}(x)$.

The EBM (associated with the GFlowNet) is trained via the Contrastive
Divergence algorithm [4], where the negative samples $x^{\prime}$ are sampled
from a GFlowNet policy ($x^{\prime}\sim P_{T}(x^{\prime})$). Note that,
$P_{T}(x^{\prime})$ denotes marginalizing the forward trajectory distribution
$P_{F}(\tau)$ over its non-terminating states. With contrastive divergence,
one samples $x^{\prime}$ from a K-step Metropolis Hastings chain, which starts
from a true data sample $x\sim p_{data}(x)$. Thanks to the backward policy
$P_{B}$, one can easily perform such MCMC steps with a GFlowNet [29]. For a
perfectly trained GFlowNet, the terminating probability distribution will be
equal to the energy distribution. The updates to $\phi$ (the parameter of the
energy model) are made proportional to Eq. 2:

$\mathbb{E}_{x\sim
p_{data}(x)}\nabla_{\phi}\mathcal{E}_{\phi}(x)-\mathbb{E}_{x^{\prime}\sim
P_{T}(x^{\prime})}\nabla_{\phi}\mathcal{E}_{\phi}(x^{\prime}).$ (2)

We introduce JEB-GFN, which extends the above formulation to model the joint
distribution over multiple variables, in particular, data $x$ and
corresponding label $y$. The unnormalized joint probability distribution is
denoted by $\mathcal{E}_{\phi}(x,y)$. The parameters of the EBM are updated
proportional to

$\mathbb{E}_{(x,y)\sim
p_{data}(x,y)}\nabla_{\phi}\mathcal{E}_{\phi}(x,y)-\mathbb{E}_{(x^{\prime},y^{\prime})\sim
P_{T}(x^{\prime},y^{\prime})}\nabla_{\phi}\mathcal{E}_{\phi}(x^{\prime},y^{\prime})$
(3)

Similar to the EB-GFN formulation, the EBM is trained via contrastive
divergence, where a K-Step MCMC is used to generate negative samples, given
data sampled from the original data distribution. This new generated sample is
accepted or rejected depending on the Metropolis-Hastings rule.

Concretely, first, a K-Step backward trajectory ($\tau$) is sampled with a
multiple calls to the backward policy $P_{B}$, from $(x,y)$ to create an
intermediate data sample / state say $(x_{int},y_{int})$. From this
intermediate data sample, a forward trajectory is sampled ($\tau^{\prime}$)
via a forward policy $P_{F}$, to a terminal state $(x^{\prime},y^{\prime})$.
Note that the transition probability from state $x$ to state $x^{\prime}$ is
defined as $P_{B}(\tau|x,y)P_{F}(\tau^{\prime})$, while the reverse transition
probability is defined as
$P_{B}(\tau^{\prime}|x^{\prime},y^{\prime})P_{F}(\tau)$. The transition of
$(x,y)$ to $(x^{\prime},y^{\prime})$ is accepted using the Metropolis-Hastings
rule, with the probability

$A_{\tau,\tau^{\prime}}((x,y)\xrightarrow[]{}(x^{\prime},y^{\prime}))\overset{\Delta}{=}min[1,\frac{e^{-\mathcal{E}_{\phi}(x^{\prime},y^{\prime})}P_{B}(\tau|x,y)P_{F}(\tau^{\prime})}{e^{-\mathcal{E}_{\phi}(x,y)}P_{B}(\tau^{\prime}|x^{\prime},y^{\prime})P_{F}(\tau)}].$
(4)

Note that the EBM and the GFlowNet are trained via an alternating optimization
routine. The procedure for joint training of the GFlowNet and EBM is described
in Algorithm 1.

Input: $D=\\{x_{i},y_{i}\\},i=1,\dots,N$: Dataset (size $N$) consisting of
data: $x_{i}$ corresponding to labels: $y_{i}$

Initialization: Initalize GFlowNet’s $P_{F}$, $P_{B}$ & $Z$ with parameters
$\theta$ and the energy function $\mathcal{E}_{\phi}$ with parameters $\phi$

for  _$i=1$ to $T$_ do

$\bullet$ Sample forward trajectory $\tau\sim P_{F}(\tau)$

$\bullet$ Compute $P_{B}(\tau)$ given $\tau$, where $\tau$ is generated from
the previous step

$\bullet$ Update $\theta$ via the gradient step on
$\mathcal{L}_{\theta}(\tau)$ (Eq. 1) with the reward function being
$e^{-\mathcal{E}_{\phi}(x,y)}$

$\bullet$ Generate horizon $K$ depending on iteration $i$

$\bullet$ Sample a $K$-step backward trajectory from
$P_{B}(\cdot|\cdot;\theta)$: $\tau$ = ($[x,y]$ = $S_{D}$ $\rightarrow$
$S_{D-1}$ $\rightarrow$ $...$ $\rightarrow$ $S_{D-K})$

$\bullet$ Sample a $K$-step forward trajectory from
$P_{F}(\cdot|\cdot;\theta)$: $\tau$ = ($S_{D-K}$ $\rightarrow$
$S_{D-K+1}^{\prime}$ $\rightarrow$ $...$ $\rightarrow$ $S_{D}^{\prime}$ $=$
$[x^{\prime},y^{\prime}]$)

$\bullet$ Accept or reject $[x^{\prime},y^{\prime}]$ via Eq. 4; set
$[x^{\prime},y^{\prime}]$ = $[x,y]$ if reject.

$\bullet$ Update $\phi$ with gradient
$\mathcal{E}_{\phi}(x,y)-\mathcal{E}_{\phi}(x^{\prime},y^{\prime})$

end for

Algorithm 1 JEB-GFN joint training framework

##  4 Empirical Results

###  4.1 Synthetic Data

In this section, we first analyse JEB-GFNs on a set of synthetic tasks adopted
from [29]. To motivate the discussion, we first demonstrate a failure mode of
naive conditional generation with GFlowNet. We then demonstrate that JEB-GFN
learns the joint $p(x,y)$ and can generate accurate conditional samples
$p(x|y)$.

Pitfalls of the Naive Conditional Generation Baseline

One way to sample from the joint $p(x,y)$ would be to train GFlowNets
conditioned on $y$, with the reward being the probability output by a
classifier $p(y|x)$, trained on the original dataset $\mathcal{D}$. We
evaluate this approach, called CGFN, on the two moons
dataset111https://scikit-
learn.org/stable/modules/generated/sklearn.datasets.make_moons.html. The
procedure is as follows: train a binary classifier on the two moons dataset.
Next we train two GFlowNets, one for each class, to generate samples belonging
to a particular class. The reward for each of these GFlowNets, was the
confidence of the classifier for the generated data for the corresponding
class. Following [29] we convert these 2D floating point values to 16-bit
grey-code. The GFlowNets and the Classifier were trained on this 16-bit grey
code data. Hence each of these GFlowNets generates a 32 bit sequence
autoregressively. At each step the GFN outputs whether a bit at a particular
location is either 1 or 0. The samples generated by the GFlowNets are
visualized in Figure 1(a) and Figure 1(b). As one can observe, the samples
generated by such an approach are very noisy and do not resemble the original
data distribution. A hypothesis to explain this performance is that neural
network classifiers tend to be overconfident on out-of-sample data.
Consequently, this results in the GFlowNet modelling the "wrong" reward.

We then trained a JEBGFN on the same setting (i.e. the two moons dataset).
Again, these 2D floating point values were converted to 16-bit grey-code. In
this case since the JEBGFN outputs both the data, and the class it belongs to,
the JEBGFN generates a 33 bit sequence (32 bits for the data and one bit for
the class). Again this sequence is generated autoregressively. The samples
generated by the GFN and the learned conditional unormalised densities
(conditional EBMs) are visualized in Figure 1(c) and Figure 1(d) respectively.
Clearly, the JEBGFN is able to learn the original data distribution, given
data samples. In sharp contrast to naive approach described above, as JEB-GFN
estimates the unnormalised probability density via contrastive divergence,
which results in a low density in out-of-sample regions as expected.

(a) CGFN samples conditioned on label 0

(b) CGFN samples conditioned on label 1

(c) Samples generated by JEB-GFN. Red and blue denote conditioning on labels 0
and 1 respectively

(d) Energy Function learned by JEB-GFN

Figure 1: (Top) shows a GFlowNet trained with a classifier as reward (CGFN),
cannot estimate the underlying distribution. (Bottom) shows that a JEB-GFN can
accurately model the underlying data distribution.

####  4.1.1 Four Gaussians

We now consider a synthetic multi-class task where the data is independently
sampled from four 2D Gaussian Distributions, and was assigned a different
label depending on the index of the Gaussian that the data was sampled from.
The four 2D Gaussians were characterised by the following means:
$(2\sqrt{2},0),(-2\sqrt{2},0),(0,2\sqrt{2}),(0,-2\sqrt{2})$ and each had a
standard deviation of $0.5$. Similar to the previous experiment, the sampled
2D floating point values were converted to 16-bit grey-code. The JEBGFN
autoregressively generates a 34 bit sequence (32 bits for the data, 2 bits for
the class). The results can be observed in Figure 2. We can observe that the
JEBGFN is indeed able to conditionally sample data given a particular class.
Moreover, the conditional unnormalised density shows that the underlying
distribution is being learned faithfully.

(a) Generated Samples

(b) Energy Function learned with JEB-GFN

Figure 2: (Left) depicts the samples generated from JEB-GFN, conditioned on
each label. Each color represents a different label. (Right) The conditional
energy function (conditioned on the label).

####  4.1.2 MNIST

The prior two experiments show the efficacy of our method on a dataset where
the probability density associated with the dataset is already supplied
apriori, i.e. we have the ability to potentially sample infinite datapoints
from the data distribution. Two questions that remain unanswered are a) How
does our method perform in a limited data setting b) How does our method
perform on higher dimensional datasets. Through this experiment we try to
answer both these questions. The dataset created for the experiment was
obtained by sampling the following digits from the MNIST dataset [9]:
$0,1,2,3$. Note that since we don’t have the "true" probability density from
where the images were sampled from we work only with a limited set of samples
(the data points provided in the original dataset). Following [29] we binarize
the original MNIST Data (converting each pixel value to either a 1 or a 0),
and then convert it to a 1D array (by concatenating all the rows). Similar to
the prior experiments, even in this case the JEB-GFN autoregressively
generates a bit sequence, however in this case the length is 786 bits (784
bits for the data and 2 bits for the class) The results of the samples
obtained from the GFlowNet can be observed in Figure 3. One can observe that
the JEBGFN is able to generate samples close to the original data
distribution. Note that, these results will not be exactly the same as the
original dataset because the data has been binarized.

(a) Samples conditioned on label $0$

(b) Samples conditioned on label $1$

(c) Samples conditioned on label $2$

(d) Samples conditioned on label $3$

Figure 3: Samples generated with JEB-GFN trained on MNIST

###  4.2 Generating Antimicrobial Peptides

####  4.2.1 Single Round

The goal of this experiment was to generate high scoring peptides (short
protein sequences) where the scores are determined based on their anti-
microbial properties. We use the dataset from [12], with 3219 AMPs and 4611
non-AMP sequences obtained from the DBAASP database [22]. Note that the
sequences considered were of length lower than $50$, and the vocabulary size
considered was of size 20 (20 amino acids). The JEBGFN autoregressively
generates a tuple $(x,y)$, having a maximum combined length of 51. (50 for $x$
and 1 for the class $y$ which indicates whether a give sequence is
antimicrobial). All sequences having a score greater than $0.75$ are assigned
the label $1$ corresponding to the positive class, and all the others are
assigned the label corresponding to the negative class $0$. The sequences are
evaluated by an oracle. Note that, in this case, this is not a bit sequence
but a sequence of integers. The range of the integers lie between 0-19,
(representing 20 amino acids). Each of these integers are converted to a 256
dimensional learnt embedding, where the position information is also
incorporated into this embedding space. We use a 3-layer MLP with 256 units in
each hidden layer for both the GFLowNet and EBM. We use the Adam optimizer to
train the JEBGFN. We search the learning rate of the GFlowNet in
$\\{10^{-4},10^{-3},10^{-2}\\}$, the learning rate of the energy function in
$\\{10^{-4},10^{-3},10^{-2}\\}$. The coefficient of ${l}_{2}$ regularization
for both the GFlowNet and the EBGFN was searched in
$\\{0,0.01,0.02,0.05,0.1\\}$. We compare our method against the following
baseline methods: a) GFlowNets [2] b) GDGDS [18]: a goal-directed conditional
generative model , c) SSVAE [13]: a conditional molecule design method, based
on a learned semi-supervised VAE. We use the performance (average reward of
Top $K$ sequences), diversity (average pairwise distance among the Top $K$
sequences) and novelty (average pairwise distance between Top $K$ generated
sequences and the initial dataset) metrics from [12].

We can observe that our model performs significantly better than the other
baselines across all the metrics. A major point to notice is that our method
generates sequences that are more diverse than all the other baselines. Note
that this experiment took $\sim 0.06$ GPU Day(s) to complete, on NVIDIA A100
GPUs.

Table 1: Single round results on the AMP Task with $K=100$.  | Performance | Diversity | Novelty
---|---|---|---
JEB-GFN | $\mathbf{0.858\pm 0.022}$ | $\mathbf{33.50\pm 0.25}$ | $\mathbf{24.25\pm 0.90}$
GFlowNet | $0.753\pm 0.023$ | $16.43\pm 4.11$ | $18.38\pm 3.09$
GDGDS | $0.712\pm 0.009$ | $11.56\pm 2.71$ | $12.34\pm 4.12$
SSVAE | $0.743\pm 0.012$ | $14.75\pm 1.98$ | $11.31\pm 1.31$

####  4.2.2 Active Learning

Through this experiment we try to understand the efficacy of JEB-GFN in the
active learning setting, for the AMP Generation task from [12]. We sample
$1000$ anti microbial peptides after each round of active learning and
evaluate the scores by an oracle. All sequences having a score greater than
$0.75$ are assigned the label $1$ corresponding to the positive class, and all
the others are assigned the label corresponding to the negative class $0$.
Note that this continues for 10 rounds of active learning. After 10 rounds of
active learning we pick the Top $K$ sequences (here $K=100$) and then evaluate
them. We use the same model details and hyperparameters from Section 4.2.1. We
compare our method against the following baseline methods: a) GFlowNet-AL b)
DynaPPO [1]: A method that uses model-based RL algorithm, for designing
biological sequences. c) COMs [25]: a method that relies on generating novel
candidates by optimizing known candidates against a learned conservative
model. d) GDGDS e) SSVAE. We use the metrics from [12], as in Section 4.2.1.
The results for the same can be found in Table 2. While the highest score
(performance) is obtained by the DynaPPO method, it severly underperforms on
the diversity and novelty metrics. Moreoever, though JEB-GFN performs worse
than GFlowNet-AL and DynaPPO (in terms of the performance), it does
significantly better than all other methods in terms of its diversity ($\sim
1.5$x better than its closest baseline). Considering the diverse nature of
viable biological targets, we believe that a method producing sequences having
a much higher diversity, while also having a high score, would have a higher
utility in practice. Note that this experiment took $\sim 1$ GPU Day(s) to
complete, on NVIDIA A100 GPUs.

Table 2: Active Learning Results on the AMP Task with $K=100$.  | Performance | Diversity | Novelty
---|---|---|---
JEB-GFN | $0.922\pm 0.004$ | $\mathbf{31.83\pm 0.837}$ | $24.26\pm 0.965$
GFlowNet-AL | $\mathbf{0.932\pm 0.002}$ | $22.34\pm 1.24$ | $\mathbf{28.44\pm 1.32}$
DynaPPO | $\mathbf{0.938\pm 0.009}$ | $12.12\pm 1.71$ | $9.31\pm 0.69$
COMs | $0.761\pm 0.009$ | $19.38\pm 0.14$ | $26.47\pm 1.3$
GDGDS | $0.785\pm 0.021$ | $9.41\pm 1.72$ | $10.73\pm 2.11$
SSVAE | $0.797\pm 0.009$ | $10.47\pm 1.32$ | $10.77\pm 1.39$

##  5 Conclusion and Future Work

Motivated from fact that generating viable biological sequences is a time
consuming and laborious process, we extend EB-GFN to learn the joint over
multiple variables. We validated the efficacy of our method over both
synthetic and real-world tasks. Moreover, we have also shown the effectiveness
of JEB-GFN in the active learning setting. One limitation of the above method
is that the sequences generated can only consist of discrete values. We
believe that modelling the joint over both discrete and continuous values
would lead to the generation of sequences (conditionally) with an increased
level of granularity, which we leave as future work. We do not foresee any
negative societal impact from this work but as with all machine learning
algorithms it could potentially be used by nefarious agents.

## References

  * [1] C. Angermueller, D. Dohan, D. Belanger, R. Deshpande, K. Murphy, and L. Colwell.  Model-based reinforcement learning for biological sequence design.  In International conference on learning representations, 2019.
  * [2] E. Bengio, M. Jain, M. Korablyov, D. Precup, and Y. Bengio.  Flow network based generative models for non-iterative diverse candidate generation.  Advances in Neural Information Processing Systems, 34, 2021.
  * [3] Y. Bengio, T. Deleu, E. J. Hu, S. Lahlou, M. Tiwari, and E. Bengio.  Gflownet foundations.  arXiv preprint arXiv: Arxiv-2111.09266, 2021.
  * [4] M. A. Carreira-Perpinan and G. Hinton.  On contrastive divergence learning.  In International workshop on artificial intelligence and statistics, pages 33–40. PMLR, 2005.
  * [5] T. Che, R. Zhang, J. Sohl-Dickstein, H. Larochelle, L. Paull, Y. Cao, and Y. Bengio.  Your gan is secretly an energy-based model and you should use discriminator driven latent sampling.  Advances in Neural Information Processing Systems, 33:12275–12287, 2020.
  * [6] Z. Chen, Z. Yang, X. Wang, X. Liang, X. Yan, G. Li, and L. Lin.  Multivariate-information adversarial ensemble for scalable joint distribution matching.  In International Conference on Machine Learning, pages 1112–1121. PMLR, 2019.
  * [7] V. Chenthamarakshan, P. Das, S. C. Hoffman, H. Strobelt, I. Padhi, K. W. Lim, B. Hoover, M. Manica, J. Born, T. Laino, and A. Mojsilovic.  Cogmol: Target-specific and selective drug design for covid-19 using deep generative models.  arXiv: Learning, 2020.
  * [8] P. Das, T. Sercu, K. Wadhawan, I. Padhi, S. Gehrmann, F. S. Cipcigan, V. Chenthamarakshan, H. Strobelt, C. D. Santos, P.-Y. Chen, Y. Y. Yang, J. P. K. Tan, J. L. Hedrick, J. Crain, and A. Mojsilovic.  Accelerated antimicrobial discovery via deep generative models and molecular dynamics simulations.  Nature Biomedical Engineering, 5:613–623, 2021.
  * [9] L. Deng.  The mnist database of handwritten digit images for machine learning research.  IEEE Signal Processing Magazine, 29(6):141–142, 2012.
  * [10] Y. Du and I. Mordatch.  Implicit generation and generalization in energy-based models.  arXiv preprint arXiv:1903.08689, 2019.
  * [11] T. He, B. McCann, C. Xiong, and E. Hosseini-Asl.  Joint energy-based model training for better calibrated natural language understanding models.  In Proceedings of the 16th Conference of the European Chapter of the Association for Computational Linguistics: Main Volume, pages 1754–1761, 2021.
  * [12] M. Jain, E. Bengio, A.-H. Garcia, J. Rector-Brooks, B. F. P. Dossou, C. Ekbote, J. Fu, T. Zhang, M. Kilgour, D. Zhang, L. Simine, P. Das, and Y. Bengio.  Biological sequence design with gflownets.  arXiv preprint arXiv: Arxiv-2203.04115, 2022.
  * [13] S. Kang and K. Cho.  Conditional molecular design with deep generative models.  Journal of chemical information and modeling, 59(1):43–52, 2018\.
  * [14] J. Kelly, R. Zemel, and W. Grathwohl.  Directly training joint energy-based models for conditional synthesis and calibrated prediction of multi-attribute data.  arXiv preprint arXiv:2108.04227, 2021.
  * [15] Y. LeCun, S. Chopra, R. Hadsell, M. Ranzato, and F. Huang.  A tutorial on energy-based learning.  Predicting structured data, 1(0), 2006.
  * [16] S. Li, Y. Du, G. M. van de Ven, and I. Mordatch.  Energy-based models for continual learning.  arXiv preprint arXiv:2011.12216, 2020.
  * [17] N. Malkin, M. Jain, E. Bengio, C. Sun, and Y. Bengio.  Trajectory balance: Improved credit assignment in gflownets.  arXiv preprint arXiv: Arxiv-2201.13259, 2022.
  * [18] A. Mollaysa, B. Paige, and A. Kalousis.  Goal-directed generation of discrete structures with conditional generative models.  Advances in Neural Information Processing Systems, 33:21923–21933, 2020.
  * [19] H. Moss, D. Leslie, D. Beck, J. Gonzalez, and P. Rayson.  Boss: Bayesian optimization over string spaces.  Advances in neural information processing systems, 33:15476–15486, 2020.
  * [20] M. M. Mullis, I. M. Rambo, B. J. Baker, and B. K. Reese.  Diversity, ecology, and prevalence of antimicrobials in nature.  Frontiers in microbiology, 10:2518, 2019.
  * [21] W. Nie, A. Vahdat, and A. Anandkumar.  Controllable and compositional generation with latent-space energy-based models.  Advances in Neural Information Processing Systems, 34, 2021.
  * [22] M. Pirtskhalava, A. A. Amstrong, M. Grigolava, M. Chubinidze, E. Alimbarashvili, B. Vishnepolsky, A. Gabrielian, A. Rosenthal, D. E. Hurt, and M. Tartakovsky.  Dbaasp v3: database of antimicrobial/cytotoxic activity and structure of peptides as a resource for development of new therapeutics.  Nucleic acids research, 49(D1):D288–D297, 2021.
  * [23] Y. Pu, S. Dai, Z. Gan, W. Wang, G. Wang, Y. Zhang, R. Henao, and L. C. Duke.  Jointgan: Multi-domain joint distribution learning with generative adversarial nets.  In International Conference on Machine Learning, pages 4151–4160. PMLR, 2018.
  * [24] Y. Song and D. P. Kingma.  How to train your energy-based models.  arXiv preprint arXiv:2101.03288, 2021.
  * [25] B. Trabucco, X. Geng, A. Kumar, and S. Levine.  Design-bench: Benchmarks for data-driven offline model-based optimization.  arXiv preprint arXiv:2202.08450, 2022.
  * [26] C. M. Van Oort, J. B. Ferrell, J. M. Remington, S. Wshah, and J. Li.  Ampgan v2: Machine learning-guided design of antimicrobial peptides.  Journal of chemical information and modeling, 61(5):2198–2207, 2021\.
  * [27] J. Wang, C.-Y. Hsieh, M. Wang, X. Wang, Z. Wu, D. Jiang, B. Liao, X. Zhang, B. Yang, Q. He, D. Cao, X. Chen, and T. Hou.  Multi-constraint molecular generation based on conditional transformer, knowledge distillation and reinforcement learning.  Nat. Mach. Intell., 3:914–922, 2021.
  * [28] B. Xie, L. Yuan, S. Li, C. H. Liu, X. Cheng, and G. Wang.  Active learning for domain adaptation: An energy-based approach.  arXiv preprint arXiv:2112.01406, 2021.
  * [29] D. Zhang, N. Malkin, Z. Liu, A. Volokhova, A. Courville, and Y. Bengio.  Generative flow networks for discrete probabilistic modeling.  arXiv preprint arXiv: Arxiv-2202.01361, 2022.


\end{bmatrix}\right).
\end{equation}

Noting that (<ref>) is separable along the components of $y$ in $\R^{E_1}$ and the components of $y$ in $\R^{E_2}$ shows that it is equivalent to (<ref>).

For the second diagram, note that the node incidence matrix of an empty flow network is necessarily the unique linear map $\R^0\to\R^0$. Similarly, for an empty flow network $b$ must be $0\in\R^0$. With this case in mind, following the top and bottom paths both yield the constant zero function, as desired.

With this algebra morphism established, we can finally turn to the composite given by the following diagram:

\begin{equation}\label{eq:big_daddy}\begin{tikzcd}
\FG && \Conc && \NDD && {\NDD_D}
\arrow["{\netflow}", Rightarrow, from=1-1, to=1-3]
\arrow["{\superg}", Rightarrow, from=1-3, to=1-5]
\arrow["{\Euler^\gamma}", Rightarrow, from=1-5, to=1-7]
\end{tikzcd}\end{equation}

This transformation takes a flow network $G$ to the discrete non-deterministic dynamical system that solves the minimum cost network flow problem on $G$ using supergradient ascent on the dual problem $q_G$. The results proven throughout this paper show that this composite respects any decomposition of $G$. In particular, if we completely decompose $G$ with a UWD where each inner box contains only a single edge, then applying the composite in (<ref>) yields the standard dual decomposition algorithm as in Section 2.A of [48].

This decomposition algorithm is flat, in that it decomposes into a single master problem and scalar subproblems for every eadge. More interestingly, we can apply (<ref>) to other decomposition structures of $G$. For example, we could decompose $G$ by densely connected components. Each connected component can then be fully decomposed by edges to yield a hierarchical decomposition. The subsequent section illustrates how exploiting such higher-level compositional structure can be beneficial.

§ IMPLEMENTATION AND NUMERICAL RESULTS

This section discusses our prototype implementation of this framework in a Julia package called AlgebraicOptimization.jl[<https://github.com/AlgebraicJulia/AlgebraicOptimization.jl>]. We used this package to perform numerical experiments demonstrating that the exploitation of higher-level compositional structure of flow networks can result in a faster solver than standard dual decomposition for MCNF.

§.§ The AlgebraicOptimization Software Architecture

The framework developed in this paper consists of a collection of gradient-based algebra morphisms

\begin{equation}
g\maps O\Rightarrow D,
\end{equation}

where $O$ is one of the optimization problem algebras and $D$ is one of the dynamical system algebras (visualized in Figure <ref>). These algebra morphisms take optimization problems as inputs and produce iterative solution algorithms as outputs. Thus for each $g\maps O\Rightarrow D$, we require the following:

* an implementation of $O$ as a UWD-algebra,

* an implementation of $D$ as a UWD-algebra,

* an implementation of the components of $g$.

The Catlab.jl software package for categorical computing [49] provides a generic interface for specifying UWD-algebras, while AlgebraicDynamics.jl implements the $\Dynam$ and $\Dynam_D$ UWD-algebras as well as the $\Euler^\gamma$ algebra morphism [17]. So, for example, to specify $\gd^\gamma\maps \Opt\Rightarrow \Dynam_D$, we just needed to specify $\Opt$ as a UWD-algebra using Catlab and the components of $\flow\maps\Opt\Rightarrow\Dynam$, which can then be composed with $\Euler^\gamma$ to recover $\gd^\gamma$. The gradient flow dynamical system for a given objective function is straightforward to implement thanks to Julia's extensive support for automatic differentiation with packages such as ForwardDiff.jl [50].

The ease with which we can implement this framework illustrates one of the benefits of formalizing these optimization problems with category theory: with a package such as Catlab for categorical computing, the translation of mathematical definitions to usable software is straightforward. Moreover, thanks to the intuitive graphical synatx of UWDs, an application-oriented user of AlgebraicOptimization need not be concerned with the underlying category theory
to benefit from these results. They can instead simply specify a UWD and subproblems to fill each inner box, and AlgebraicOptimization will automatically generate a
decomposed solution method for their problem. The user can also specify whether the solution method is serial (distribute to subsystems, update subsystems sequentially, then collect from subsystems), or parallel (distribute to subsystems, update subsystems in parallel, then collect from subsystems).

§.§ Compositional Network Flow Experiments

Benchmark results of standard dual decomposition and hierarchical dual decomposition on random composite flow networks. The left column reports results for varying the number of nodes in each components subgraph from 10 to 150 while holding connectivity fixed at 0.2. The right column reports results for varying the graph connectivity parameter of each subgraph from 0.1 to 1.0 while holding the number of nodes per subgraph fixed at 80. The top row compares execution times, the second row compares total memory used, and the bottom row charts the ratios between standard and hierarchical dual decomposition for each of these quantities.

Intuitively, we expect that compositional modeling of MCNF problems will be practically useful because many real world applications of network flow involve flow networks with several densely connected components and sparse connectivity among these components. For example, a traffic flow network may model several densely connected urban areas coupled by relatively few interstates and highways. This is an example of hierarchical compositional structure where the first level decomposes each urban area from the others and the second level decomposes the edges within each urban area.

Motivated by such applications, we conducted a brief experimental study to measure the performance of standard dual decomposition against the hierarchical dual decomposition enabled by our framework when solving the minimum cost network flow problem. We worked with flow networks whose compositional structure was given by the UWD $\Phi\maps [3]+[2]+[2]\to\emptyset$ in Figure <ref>. The graphs were randomly generated to avoid implicitly favoring any particular solution method. More specifically, to generate test data for this algorithm, we built three component flow networks $(G_1,G_2,G_3)$ by starting with three connected Erdős-Rényi random graphs [51], assigning to each edge a random convex quadratic cost function, and to the vertices a random vector of inflows and outflows which summed to 0. We then randomly designated 3 boundary vertices for $G_1$ and 2 boundary vertices for each of $G_2$ and $G_3$. Thus constructed, $G_1$ fills the box of $\Phi$ with three ports, and $G_2$ and $G_3$ fill the boxes with two ports. Given $(G_1, G_2, G_3)$, we examined two ways to proceed towards solving the corresponding MCNF problem:

* (Standard Dual Decomposition) This algorithm is given by first composing component subgraphs with $\FG$ and then taking the dual decomposition optimizer of the composite graph.

* (Hierarchical Dual Decomposition) This algorithm is given by taking the dual decomposition optimizer of each component subgraph and composing these systems with $\Dynam_D$.

Standard dual decomposition corresponds to following the top path of the following diagram while hierarchical dual decomposition corresponds to following the bottom path.

\[\begin{tikzcd}
{\FG([3])\times\FG([2])\times\FG([2])} &&& {\FG(\emptyset)} \\
{\DiffConc([3])\times\DiffConc([2])\times\DiffConc([2])} &&& {\DiffConc(\emptyset)} \\
{\Dynam([3])\times\Dynam([2])\times\Dynam([2])} &&& {\Dynam(\emptyset)} \\
{\Dynam_D([3])\times\Dynam_D([2])\times\Dynam_D([2])} &&& {\Dynam_D(\emptyset)}
\arrow["{\FG(\Phi)}", from=1-1, to=1-4]
\arrow["{\netflow_\emptyset}", from=1-4, to=2-4]
\arrow["\flow_\emptyset", from=2-4, to=3-4]
\arrow["{\Euler^\gamma_\emptyset}", from=3-4, to=4-4]
\arrow["{\netflow_{[3]}\times\netflow_{[2]}\times\netflow_{[2]}}"', from=1-1, to=2-1]
\arrow["{\flow_{[3]}\times\flow_{[2]}\times\flow_{[2]}}"', from=2-1, to=3-1]
\arrow["{\Euler^\gamma_{[3]}\times\Euler^\gamma_{[2]}\times\Euler^\gamma_{[2]}}"', from=3-1, to=4-1]
\arrow["{\Dynam_D(\Phi)}"', from=4-1, to=4-4]
\end{tikzcd}\]

Note that for all experiments performed, we used strictly convex cost functions for each edge, implying that the dual function of a flow network generated by the $\netflow$ transformation is differentiable. We were therefore able to use standard gradient flow instead of supergradient flow.
Also recall that the algebras in this diagram represent the UWD-algebras of the corresponding finset algebras defined in this paper. So for example, $\FG(\emptyset)$ represents the set of all open flow networks with empty boundary, rather than the empty flow network.
Since the above diagram is the naturality square of a composite of natural transformations, the dynamical systems generated by following the top and bottom paths are equivalent in that they follow the same trajectory.

Given these two methods of building solvers from problem data, a fair comparison between them must also account for relevant properties of the random graphs. In particular, one might expect that the size, i.e. number of nodes and edges, in the component graphs would affect execution time and memory usage of each algorithm. As such we sought to compare both standard and hierarchical dual decomposition using graphs of various sizes. More specifically, we conducted two experiments to compare standard dual decomposition and hierarchical dual decomposition. First note that an Erdős-Rényi random graph $G(N,p)$ is parameterized by two values: $N\geq 0$ is the number of vertices in the graph and $0\leq p\leq1$ represents the probability that a given edge is included in the graph. Therefore, increasing $p$ increases the connectivity of $G(N,p)$. For the first experiment, we held $p$ at the fixed value of 0.2 and varied the number of nodes for each component subgraph. Conversely, for the second experiment, we held the number of nodes of each component subgraph fixed at 80 and varied $p$. For each selection of parameters, we measured the elapsed time and total memory consumption to run 10 outer iterations of dual decomposition. We used the Optim.jl package [52] to solve the inner problems of computing optimal flows for edges given fixed values of the dual variables.

The results of these experiments are summarized in Figure <ref>. These results clearly demonstrate the benefits of hierarchical dual decomposition over standard dual decomposition when hierarchical structure exists. Note that we do not compare against state of the art MCNF algorithms, but rather aim to show some practical utility of hierarchical structure.  In particular, the advantages of hierarchical decomposition are more pronounced on dense graphs with many nodes. This is expected because although both algorithms must solve the same number of scalar minimization problems (one for each edge), these minimizations can be computed more efficiently in the hierarchical case because they are only parameterized by the dual variables of the subgraph they belong to, rather than all the dual variables in the case of standard dual decomposition.

§ CONCLUSION

§.§ Summary

We have presented a novel algebraic/structural perspective on first-order optimization decomposition algorithms using the abstractions of undirected wiring diagrams (UWDs), operad algebras, and algebra morphisms.
Specifically, we demonstrated that whenever a first-order method defines an algebra morphism from a UWD-algebra of problems to a UWD-algebra of dynamical systems, that method decomposes problems defined on arbitrary UWDs and provides distributed message-passing solution algorithms. We leveraged this paradigm to
unify and generalize many known algorithms within a common framework, including gradient descent (Theorem <ref>), Uzawa's algorithm (Theorem <ref>), primal and dual decomposition (Examples <ref> and <ref>), and subgradient variants of each (Corollary <ref>). Working within this framework allowed the statement of the compositional data condition (Definition <ref>), which is a novel sufficient condition for when optimization problems parameterized by compositional data are decomposable. We subsequently showed that the minimum cost network flow problem satisfied the compositional data condition by defining an algebra morphism from a UWD-algebra of flow networks to the UWD-algebra of concave maximization problems. Experimental results using our implementation of this compositional framework highlighted the computational benefits afforded by the ability to exploit hierarchical and compositional problem structure.

§.§ Impacts

We believe this framework will have many positive impacts for both optimization theory and practice. Most immediately, our implementation of this framework provides a convenient tool for engineers to compositionally build complex optimization problems in a hierarchical fashion using the intuitive graphical syntax of undirected wiring diagrams. Such specifications are then automatically transformed into performant distributed solution algorithms by leveraging the algebra morphisms described in this paper. We also hope that optimization researchers can utilize this general unified language for optimization decomposition methods to more easily build novel algorithms for specialized applications (e.g. by finding an application for a hierarchical combination of both primal and dual decomposition). Furthermore, the compositional data condition provides a principled way of proving when problems defined by compositional data are decomposable and generating appropriate distributed solution algorithms. We hope that this condition can be used to derive hierarchical decomposition methods for many more applications beyond the network flow example presented in the paper.

In terms of theoretical impacts, this work has situated optimization decomposition methods in the broad mathematical landscape of category theory. Seen in this light, the underlying mathematical structure of optimization decomposition methods is the same as or similar to the structures underlying entropy maximization of thermostatic systems [53] and equilibrium behavior of Markov processes [16] in that both also tell stories of operadic dynamical systems. Also of interest is that backpropagation in neural networks has been shown to exhibit semantics for directed wiring diagrams [24]. In making such connections we hope to facilitate the transfer of ideas amongst optimization and other disciplines. Additionally, we believe that the category theoretic point of view will lead to the expansion of problems considered by optimization researchers. Indeed constructs already exist that allow one to generalize gradient based optimization to settings beyond differentiable functions over the real numbers [54, 35], for example to polynomials over arbitrary finite fields or boolean circuits. Our work has the potential to add composition mechanisms to instances of such non-traditional optimization problems.

§.§ Limitations and Future Work

There are many directions for future work. The first natural question is whether other, more sophisticated decomposition methods such as the alternating direction method of multipliers (ADMM) or second-order methods such as Newton's method can be expressed within our framework. We are hopeful that an affirmative answer can be given for both of these methods. For example, [8] defines a way to apply ADMM to decompose problems defined on factor graphs, which are similar to UWDs, and [48] defines a distributed Newton's method for solving minimum cost network flow.

There are also some fundamental limitations of our approach which we aim to address in future work. For example, Euler's method (and by extension gradient descent) is only an algebra morphism for a fixed step-size. However, many optimization algorithms use more sophisticated methods such as backtracking line search for computing appropriate step-sizes at each iteration. As such, we would like to develop a theory of approximate algebra morphisms and see if we can bound the approximation error based on the compositional structure of a given problem.

Another limitation of our framework is that the distributed solution algorithms generated by applying our algebra morphisms are synchronous, meaning they require a blocking synchronization step after every parallel computation step. In reality, many algorithms such as distributed gradient descent and Uzawa's can be run asynchronously while still guaranteeing convergence [55]. We thus would like to extend the existing algebraic theory of composite dynamical systems to encompass such asynchrony.

§.§ Acknowledgements

Tyler Hanks was supported by the National Science Foundation Graduate
Research Fellowship Program under Grant No. DGE-1842473. Any opinions, findings,
and conclusions or recommendations expressed in this material are those of the author(s)
and do not necessarily reflect the views of the NSF.
Matthew Hale was supported by the Air Force Office of
Scientific Research (AFOSR) under Grant No. FA9550-23-1-0120.
Matthew Klawonn was funded by a SMART SEED grant. James Fairbanks was supported by DARPA under award no. HR00112220038.

[1]
Stephen Boyd.
Distributed Optimization and Statistical Learning via the
Alternating Direction Method of Multipliers.
Foundations and Trends® in Machine Learning, 3(1):1–122,

[2]
Neal Parikh and Stephen Boyd.
Block splitting for distributed optimization.
Mathematical Programming Computation, 6(1):77–102, March 2014.

[3]
Pontus Giselsson, Minh Dang Doan, Tamás Keviczky, Bart De Schutter, and Anders
Accelerated gradient methods and dual decomposition in distributed
model predictive control.
Automatica, 49(3):829–833, March 2013.

[4]
Dinesh Krishnamoorthy, Bjarne Foss, and Sigurd Skogestad.
A primal decomposition algorithm for distributed multistage scenario
model predictive control.
Journal of Process Control, 81:162–171, 2019.

[5]
J. F. Benders.
Partitioning procedures for solving mixed-variables programming
Numer. Math., 4(1):238–252, dec 1962.

[6]
Ragheb Rahmaniani, Teodor Gabriel Crainic, Michel Gendreau, and Walter Rei.
The Benders decomposition algorithm: A literature review.
European Journal of Operational Research, 259(3):801–817, June

[7]
Stephen Boyd, Lin Xiao, Almir Mutapcic, and Jacob Mattingly.
Notes on decomposition methods.
Notes for EE364B, Stanford University, 2015.

[8]
David Hallac, Christopher Wong, Steven Diamond, Abhijit Sharang, Rok Sosiˇc,
Stephen Boyd, and Jure Leskovec.
SnapVX: A Network-Based Convex Optimization Solver.
Journal of Machine Learning Research, 18:1–5, 2017.

[9]
Jordan Jalving, Sungho Shin, and Victor M. Zavala.
A graph-based modeling abstraction for optimization: Concepts and
implementation in plasmo.jl.
Mathematical Programming Computation, 14:699 – 747, 2022.

[10]
David I. Spivak.
The operad of wiring diagrams: formalizing a graphical language for
databases, recursion, and plug-and-play circuits, May 2013.
arXiv:1305.0297 [cs, math].

[11]
John C. Baez, Kenny Courser, and Christina Vasilakopoulou.
Structured versus Decorated Cospans.
Compositionality, 4:3, September 2022.
arXiv:2101.09363 [math].

[12]
John C. Baez and Brendan Fong.
A compositional framework for passive linear networks.
Theory and Applications of Categories, 33(38):1158–1222, 2018.

[13]
John C. Baez and Blake S. Pollard.
A Compositional Framework for Reaction Networks.
Reviews in Mathematical Physics, 29(09):1750028, October 2017.
arXiv:1704.02051 [math-ph].

[14]
John C. Baez and Jason Erbele.
Categories in control.
Theory and Applications of Categories, 30(24):836–881, 2015.

[15]
Dario Stein and Richard Samuelson.
Towards a compositional framework for convex analysis (with
applications to probability theory), 2023.

[16]
John C. Baez, Brendan Fong, and Blake S. Pollard.
A compositional framework for markov processes.
Journal of Mathematical Physics, 57(3), March 2016.

[17]
Sophie Libkind, Andrew Baas, Evan Patterson, and James Fairbanks.
Operadic Modeling of Dynamical Systems: Mathematics and
Electronic Proceedings in Theoretical Computer Science,
372:192–206, November 2022.
arXiv:2105.12282 [math].

[18]
Riccardo Scattolini.
Architectures for distributed and hierarchical model predictive
control–a review.
Journal of process control, 19(5):723–731, 2009.

[19]
Mihajlo D Mesarovic, Donald Macko, and Yasuhiko Takahara.
Theory of hierarchical, multilevel, systems.
Elsevier, 2000.

[20]
Wladyslaw Findeisen, Frederic N Bailey, Mieczyslaw Brdys, Krzysztof Malinowski,
Piotr Tatjewski, and Adam Wozniak.
Control and coordination in hierarchical systems.
John Wiley & Sons, 1980.

[21]
K.J. Arrow, L. Hurwics, and H. Uzawa.
Studies in Linear and Non-linear Programming.
Standford University Press, 1958.

[22]
Yuri Nesterov.
Primal-dual subgradient methods for convex problems.
Mathematical Programming, 120:221–259, 2009.

[23]
Michael Grant, Stephen Boyd, and Yinyu Ye.
Disciplined Convex Programming.
In Leo Liberti and Nelson Maculan, editors, Global
Optimization, volume 84, pages 155–210. Kluwer Academic Publishers,
Boston, 2006.
Series Title: Nonconvex Optimization and Its Applications.

[24]
Brendan Fong, David Spivak, and Rémy Tuyéras.
Backprop as functor: A compositional perspective on supervised
In 2019 34th Annual ACM/IEEE Symposium on Logic in Computer
Science (LICS), pages 1–13. IEEE, 2019.

[25]
G. S. H. Cruttwell, Bruno Gavranović, Neil Ghani, Paul Wilson, and Fabio
Categorical foundations of gradient-based learning, 2021.

[26]
Baike She, Tyler Hanks, James Fairbanks, and Matthew Hale.
Characterizing compositionality of lqr from the categorical
Proceedings of the 62nd IEEE Conference on Decision and
Control, 2023.

[27]
Neil Ghani, Jules Hedges, Viktor Winschel, and Philipp Zahn.
Compositional game theory.
In Proceedings of the 33rd annual ACM/IEEE symposium on logic in
computer science, pages 472–481, 2018.

[28]
Robert Atkey, Bruno Gavranović, Neil Ghani, Clemens Kupke, Jérémy
Ledent, and Fredrik Nordvall Forsberg.
Compositional game theory, compositionally.
arXiv preprint arXiv:2101.12045, 2021.

[29]
Brendan Fong and David I Spivak.
Seven sketches in compositionality: An invitation to applied category
theory, 2018.

[30]
E. Riehl.
Category theory in context.
Aurora: Dover modern math originals. Dover Publications, 2017.

[31]
Steve Awodey.
Category theory.
Oxford university press, 2010.

[32]
John Baez and Mike Stay.
Physics, topology, logic and computation: a Rosetta Stone.
Springer, 2011.

[33]
R. Tyrrell Rockafellar.
Convex Analysis.
Princeton University Press, 1970.

[34]
Stephen Boyd and Lieven Vandenberghe.
Convex Optimization.
Cambridge University Press, March 2004.

[35]
Dan Shiebler.
Generalized Optimization: A First Step Towards Category
Theoretic Learning Theory.
In Pandian Vasant, Ivan Zelinka, and Gerhard-Wilhelm Weber, editors,
Intelligent Computing & Optimization, pages 525–535, Cham, 2022.
Springer International Publishing.

[36]
Brendan Fong and David I. Spivak.
Hypergraph Categories, January 2019.
arXiv:1806.08304 [cs, math].

[37]
Claudio Hermida.
Representable Multicategories.
Advances in Mathematics, 151(2):164–225, 2000.

[38]
Brendan Fong.
The Algebra of Open and Interconnected Systems, September
arXiv:1609.05382 [math].

[39]
Philipp Schiele, Eric Luxenberg, and Stephen Boyd.
Disciplined Saddle Programming, January 2023.
arXiv:2301.13427 [cs, math].

[40]
S. Boyd, J. Duchi, M. Pilanci, and L. Vandenberghe.
Notes for EE364B, Stanford University, 2022.

[41]
Mokhtar S Bazaraa, John J Jarvis, and Hanif D Sherali.
Linear programming and network flows.
John Wiley & Sons, 2011.

[42]
Lester Randolph Ford.
Flows in networks.
Princeton University Press, 2015.

[43]
Paul A. Jensen and J. Wesley Barnes.
Network flow programming.
R.E. Krieger Publishing Company, 1987.

[44]
Jeff L Kennington and Richard V Helgason.
Algorithms for network programming.
John Wiley & Sons, Inc., 1980.

[45]
James Evans.
Optimization algorithms for networks and graphs.
CRC Press, 2017.

[46]
Ravindra K Ahuja, Thomas L Magnanti, and James B Orlin.
Network flows: Theory, applications and algorithms.
Prentice-Hall, Englewood Cliffs, New Jersey, USA Arrow, KJ
(1963). Social Choice and Individual Values, Wiley, New York. Gibbard,
A.(1973).“Manipulation of Voting Schemes: A general result”,
Econometrica, 41:587–602, 1993.

[47]
David P Williamson.
Network flow algorithms.
Cambridge University Press, 2019.

[48]
Ali Jadbabaie, Asuman Ozdaglar, and Michael Zargham.
A distributed newton method for network optimization.
In Proceedings of the 48h IEEE Conference on Decision and
Control (CDC) held jointly with 2009 28th Chinese Control
Conference, pages 2736–2741, Shanghai, China, December 2009. IEEE.

[49]
Micah Halter, Evan Patterson, Andrew Baas, and James Fairbanks.
Compositional scientific computing with catlab and semanticmodels.
arXiv preprint arXiv:2005.04831, 2020.

[50]
J. Revels, M. Lubin, and T. Papamarkou.
Forward-mode automatic differentiation in Julia.
arXiv:1607.07892 [cs.MS], 2016.

[51]
Paul L. Erdős and Alfréd Rényi.
On random graphs. i.
Publicationes Mathematicae Debrecen, 1959.

[52]
Patrick Kofod Mogensen and Asbjørn Nilsen Riseth.
Optim: A mathematical optimization package for Julia.
Journal of Open Source Software, 3(24):615, 2018.

[53]
John C. Baez, Owen Lynch, and Joe Moeller.
Compositional thermostatics.
Journal of Mathematical Physics, 64(2), February 2023.

[54]
Geoffrey SH Cruttwell, Bruno Gavranović, Neil Ghani, Paul Wilson, and Fabio
Categorical foundations of gradient-based learning.
In European Symposium on Programming, pages 1–28. Springer
International Publishing Cham, 2022.

[55]
Katherine Hendrickson and Matthew Hale.
Totally asynchronous primal-dual convex optimization in blocks, 2022.

§ CATEGORY THEORY DEFINITIONS

§.§ Categories, Functors, and Natural Transformations

A category $\mathcal{C}$ consists of a collection of objects $X,Y,Z,\dots$, and a collection of morphisms (also called arrows) $f,g,h, \dots$, together with the following data:

* Each morphism has a specified domain and codomain object; the notation $f:X\rightarrow Y$ signifies that $f$ is a morphism with domain $X$ and codomain $Y$.

* Each object $X\in \mathcal{C}$ has a distinguished identity morphism $\text{id}_X:X\rightarrow X$.

* For any triple of objects $X,Y,Z\in \mathcal{C}$ and any pair of morphisms $f:X\rightarrow Y$ and $g:Y\rightarrow Z$, there is a specified composite morphism $g\circ f:X\rightarrow Z$.

This data is subject to the following two axioms:

* (unitality) For any morphism $f:X\rightarrow Y$, we have $f\circ \text{id}_X = f = \text{id}_Y\circ f$,

* (associativity) For any objects $W,X,Y,Z\in\mathcal{C}$ and morphisms $f:W\rightarrow X$, $g:X\rightarrow Y$, $h:Y\rightarrow Z$, we have $h\circ(g\circ f) = (h\circ g)\circ f$.

Given a category $\mathcal{C}$, the opposite category, denoted $\mathcal{C}^{op}$, has

* the same objects as $\mathcal{C}$,

* a morphism $f^{op}\maps Y \to X$ for every morphism $f\maps X \to Y$ in $\mathcal{C}$.

Identities are the same as those in $\mathcal{C}$ and composition is defined by composition in $\mathcal{C}$:

\begin{equation}
g^{op}\circ f^{op} = (f\circ g)^{op}.
\end{equation}

Given categories $\mathcal{C}$ and $\mathcal{D}$, a (covariant) functor $F:\mathcal{C}\rightarrow \mathcal{D}$ consists of the following maps:

* an object map $c\mapsto F(c)$ for all $c\in \textnormal{Ob}(\mathcal{C})$,

* a morphism map respecting domains and codomains, i.e., $(f: c \to c') \mapsto (F(f): F(c) \to F(c'))$ for every morphism $f\maps c\to c'$ in $\mathcal{C}$.

These maps must satisfy the following equations:

* $F(\text{id}_c) = \text{id}_{F(c)}$ (identities are preserved),

* $F(g\circ f) = F(g)\circ F(f)$ (composition is preserved),

for all objects $c$ and composable morphisms $f,g$ in $\mathcal{C}$. A contravariant functor is a functor $F\maps \CC^{op}\to \DD$ whose domain is an opposite category.

Given a pair of functors $F,G\maps \CC\to \DD$ with common domain and codomain, a natural transformation $\theta$ from $F$ to $G$, denoted $\theta\maps F\Rightarrow G$ consists of a family of morphisms in $\DD$ indexed by objects in $\CC$, namely

\begin{equation}
\theta_X\maps F(X)\to G(X)
\end{equation}

for all $X\in \CC$. These are known as the components of $\theta$. These components are required to make the following naturality squares commute for all morphisms objects $X,Y$ and morphisms $f\maps X\to Y$ in $\CC$.

\[\begin{tikzcd}
{F(X)} && {F(Y)} \\
{G(X)} && {G(Y)}
\arrow["{F(f)}", from=1-1, to=1-3]
\arrow["{G(f)}"', from=2-1, to=2-3]
\arrow["{\theta_X}"', from=1-1, to=2-1]
\arrow["{\theta_Y}", from=1-3, to=2-3]
\end{tikzcd}\]

When each component of $\theta$ is an isomorphism, we call $\theta$ a natural isomorphism.

§.§ Universal Constructions

One of the core benefits of category theory is the ability to study an object in terms of its relationships to other objects in a category. This line of reasoning leads to very general constructions that can be formulated in any category and prove useful when specialized to familiar categories.

An initial object in a category $\mathcal{C}$ is an object $i\in\mathcal{C}$ such that for all objects $c\in\mathcal{C}$, there exists a unique morphism $i\to c$ in $\mathcal{C}$.
Dually, a terminal object in $\mathcal{C}$ is an initial object in $\CC^{op}$, i.e., an object $t\in\mathcal{C}$ such that for all objects $c\in\mathcal{C}$, there exists a unique morphism $c\to t$ in $\mathcal{C}$.

* The initial object in $\Set$ is the empty set while the terminal object[It is easy to verify that initial and terminal objects (and more generally any objects defined by universal properties) are unique up to isomorphism, so we are justified in using the terminology of “the initial/terminal object"] in $\Set$ is the singleton set $\1$.

* The vector space $\R^0$ is both the initial and the terminal object in $\Vect$.

In the following definition, a span in a category simply refers to a pair of morphisms in that category with common domain. Dually, a cospan in a category is a pair of morphisms with shared codomain.

The product of objects $X$ and $Y$ in a category $\mathcal{C}$ is a span of the form $X\xleftarrow{\pi_X} X\times Y\xrightarrow{\pi_Y} Y$ such that for any other span of the form $X\xleftarrow{f} Z\xrightarrow{g} Y$, there exists a unique morphism $Z\xrightarrow{!} X\times Y$ such that $\pi_X\circ ! = f$ and $\pi_Y\circ ! = g$. The unique morphism is called the pairing of $f$ and $g$ and denoted by $\langle f,g\rangle$.
Dually, the coproduct of objects $X$ and $Y$ in a category $\mathcal{C}$ is the product of $X$ and $Y$ in $\mathcal{C}^{op}$. This dual nature is clearly visualized by the following pair of diagrams, where the left diagram represents a product in $\mathcal{C}$ and the right diagram represents a coproduct in $\mathcal{C}$.

\[\begin{tikzcd}
X &&&& X \\
& {X\times Y} && Z && {X+Y} && Z \\
Y &&&& Y
\arrow["{\pi_X}", from=2-2, to=1-1]
\arrow["{\pi_Y}"', from=2-2, to=3-1]
\arrow[""{name=0, anchor=center, inner sep=0}, "f"', curve={height=6pt}, from=2-4, to=1-1]
\arrow[""{name=1, anchor=center, inner sep=0}, "g", curve={height=-6pt}, from=2-4, to=3-1]
\arrow["{\exists !}"', dashed, from=2-4, to=2-2]
\arrow["{\iota_X}"', from=1-5, to=2-6]
\arrow["{\iota_Y}", from=3-5, to=2-6]
\arrow["{\exists !}", dashed, from=2-6, to=2-8]
\arrow[""{name=2, anchor=center, inner sep=0}, "f", curve={height=-6pt}, from=1-5, to=2-8]
\arrow[""{name=3, anchor=center, inner sep=0}, "g"', curve={height=6pt}, from=3-5, to=2-8]
%\arrow[shorten <=3pt, shorten >=3pt, Rightarrow, no head, from=0, to=2-2]
%\arrow[shorten <=3pt, shorten >=3pt, Rightarrow, no head, from=1, to=2-2]
%\arrow[shorten <=3pt, shorten >=3pt, Rightarrow, no head, from=2, to=2-6]
%\arrow[shorten <=3pt, shorten >=3pt, Rightarrow, no head, from=3, to=2-6]
\end{tikzcd}\]

The unique arrow $X+Y\to Z$ in the coproduct diagram is called the copairing of $f$ and $g$, and denoted $[f,g]$.

* The product in $\Set$ is the Cartesian product of sets together with the natural projection maps. The coproduct in $\Set$ is the disjoint union of sets together with the natural inclusion maps.

* The product in $\Vect$ is the Cartesian product of vector spaces while the coproduct is the direct sum of vector spaces. In finite dimensions, the product and direct sum are isomorphic, meaning $\Vect$ is a special type of category called a biproduct category.

The pullback of a cospan $X\xrightarrow{f}Z\xleftarrow{g} Y$ in a category $\CC$ is a span $X\xleftarrow{p_X} P\xrightarrow{p_Y}Y$ such that $f\circ p_X = g\circ p_Y$ and for any other span $X\xleftarrow{q_X} Q\xrightarrow{q_Y}Y$, there exists a unique morphism $Q\to P$ making the following diagram commute.

\[\begin{tikzcd}
Q \\
& P & Y \\
& X & Z
\arrow["f"', from=3-2, to=3-3]
\arrow["g", from=2-3, to=3-3]
\arrow["{p_Y}", from=2-2, to=2-3]
\arrow["{p_X}"', from=2-2, to=3-2]
\arrow["{\exists!}"', dashed, from=1-1, to=2-2]
\arrow["{q_Y}", curve={height=-12pt}, from=1-1, to=2-3]
\arrow["{q_X}"', curve={height=12pt}, from=1-1, to=3-2]
\end{tikzcd}\]

The pullback object $P$ is typically denoted $X\times_Z Y$.

Dually, the pushout of a span $X\xleftarrow{f}Z\xrightarrow{g}Y$ is a cospan $X\xrightarrow{i_X} P\xleftarrow{i_Y}Y$ such that $i_X\circ f=i_Y\circ g$ and for any other cospan $X\xrightarrow{q_X} Q\xleftarrow{q_Y} Y$, there exists a unique morphism $P\to Q$ making the following diagram commute.

\[\begin{tikzcd}
Q \\
& P & Y \\
& X & Z
\arrow["f"', tail reversed, no head, from=3-2, to=3-3]
\arrow["g", tail reversed, no head, from=2-3, to=3-3]
\arrow["{i_Y}", tail reversed, no head, from=2-2, to=2-3]
\arrow["{i_X}"', tail reversed, no head, from=2-2, to=3-2]
\arrow["{\exists!}"', dashed, tail reversed, no head, from=1-1, to=2-2]
\arrow["{q_Y}", curve={height=-12pt}, tail reversed, no head, from=1-1, to=2-3]
\arrow["{q_X}"', curve={height=12pt}, tail reversed, no head, from=1-1, to=3-2]
\end{tikzcd}\]

The pushout object $P$ is typically denoted $X+_Z Y$.

* The pullback of a cospan $X\xrightarrow{f} Z\xleftarrow{g} Y$ in $\Set$ is the set $X\times_Z Y\coloneqq \{(x,y)\in X\times Y \mid f(x)=g(y)\}$ together with the projections from $X\times Y$ restricted to $X\times_Z Y$.

* The pushout of a span $X\xleftarrow{f} Z\xrightarrow{g} Y$ in $\Set$ is the set $X+_Z Y\coloneqq X+Y/\sim$ where $\sim$ is the equivalence relation $f(z)\sim g(z)$ for all $z\in Z$. The pushout also includes the maps $i_X\maps X\to X+_Z Y$ and $i_Y\maps Y\to X+_Z Y$ making $i_X\circ f = i_Y\circ g$.

§.§ Symmetric Monoidal Categories

A symmetric monoidal category consists of a category $\CC$ equipped with a functor $\otimes\maps\CC\times \CC\to \CC$ called the monoidal product and an object $1\in\CC$ called the monoidal unit, together with the following natural isomorphisms:

* (associator) $\alpha_{X,Y,Z}\maps (X\otimes Y)\otimes Z \cong X\otimes (Y\otimes Z)$,

* (left unitor) $\lambda_X\maps 1\otimes X\cong X$,

* (right unitor) $\rho_X\maps X\otimes 1\cong X$.

This data must satisfy two families of commutative diagrams called the triangle and pentagon identities which we omit for brevity.

When a category $\CC$ has finite products and a terminal object, it automatically forms a symmetric monoidal category by taking the monoidal product of morphisms $f\maps X\to Y$ and $g\maps W\to Z$ to be $\langle f\circ \pi_X, g\circ\pi_W\rangle\maps X\times W\to Y\times Z$ as in the following diagram:

\[\begin{tikzcd}
& X && Y \\
{X\times W} && {Y\times Z} \\
& W && {Z.}
\arrow["{\pi_X}", from=2-1, to=1-2]
\arrow["{\pi_W}"', from=2-1, to=3-2]
\arrow["{\pi_Y}"', from=2-3, to=1-4]
\arrow["{\pi_Z}", from=2-3, to=3-4]
\arrow["f", from=1-2, to=1-4]
\arrow["g"', from=3-2, to=3-4]
\arrow["{\langle f\circ\pi_X,g\circ\pi_W\rangle}", dashed, from=2-1, to=2-3]
\end{tikzcd}\]

This is known as a Cartesian monoidal category. For example, $(\Set,\times, 1)$ is cartesian monoidal.

Dually, when a category has finite coproducts and an initial object, it automatically forms a symmetric monoidal category by taking the monoidal product of morphisms $f\maps X\to Y$ and $g\maps W\to Z$ to be $[\iota_Y\circ f,\iota_Z\circ g]$ as in the following diagram:

\[\begin{tikzcd}
X && Y \\
& {X+W} && {Y+Z.} \\
W && {Z}
\arrow["{\iota_X}"', from=1-1, to=2-2]
\arrow["{\iota_W}", from=3-1, to=2-2]
\arrow["{[\iota_Y\circ f,\iota_Z\circ g]}", dashed, from=2-2, to=2-4]
\arrow["f", from=1-1, to=1-3]
\arrow["g"', from=3-1, to=3-3]
\arrow["{\iota_Y}", from=1-3, to=2-4]
\arrow["{\iota_Z}"', from=3-3, to=2-4]
\end{tikzcd}\]

This is known as a co-Cartesian monoidal category. For example, $(\FinSet, +, \emptyset)$ is co-Cartesian monoidal.

Given symmetric monoidal categories $(\CC, \otimes, 1)$ and $(\DD, \boxtimes, e)$, a lax symmetric monoidal functor $(\CC, \otimes, 1)\to (\DD,\boxtimes, e)$ consists of a functor $F\maps \CC\to \DD$ together with the following comparison maps:

* a morphism $\varphi_0\maps e\to F(1)$ called the unit comparison,

* a natural transformation called the product comparison with components

\begin{equation}
\varphi_{X,Y}\maps F(X)\boxtimes F(Y)\to F(X\otimes Y)
\end{equation}

for all $X,Y\in \CC$.

This data must satisfy the following axioms.

* (associativity) For all objects $X,Y,Z\in\CC$, the following diagram commutes.

\[\begin{tikzcd}
{(F(X)\boxtimes F(Y))\boxtimes F(Z)} && {F(X)\boxtimes (F(Y)\boxtimes F(Z))} \\
{F(X\otimes Y)\boxtimes F(Z)} && {F(X)\boxtimes F(Y\otimes Z)} \\
{F((X\otimes Y)\otimes Z)} && {F(X\otimes (Y\otimes Z))}
\arrow["{\alpha^\DD_{FX,FY,FZ}}", from=1-1, to=1-3]
\arrow["{\varphi_{X,Y}\boxtimes\id}"', from=1-1, to=2-1]
\arrow["{\varphi_{X\otimes Y, Z}}"', from=2-1, to=3-1]
\arrow["{F(\alpha^\CC_{X,Y,Z})}"', from=3-1, to=3-3]
\arrow["{\id\otimes\varphi_{Y,Z}}", from=1-3, to=2-3]
\arrow["{\varphi_{X,Y\otimes Z}}", from=2-3, to=3-3]
\end{tikzcd}\]

* (unitality) For all $X\in \CC$, the following pair of diagrams commutes.

\[\begin{tikzcd}
{e\boxtimes F(X)} && {F(1)\boxtimes F(X)} && {F(X)\boxtimes e} && {F(X)\boxtimes F(1)} \\
{F(X)} && {F(1\otimes X)} && {F(X)} && {F(X\otimes 1)}
\arrow["{\lambda_{F(X)}^\DD}"', from=1-1, to=2-1]
\arrow["{\varphi_0\boxtimes\id}", from=1-1, to=1-3]
\arrow["{\varphi_{1,X}}", from=1-3, to=2-3]
\arrow["{F(\lambda_X^\CC)}", from=2-3, to=2-1]
\arrow["{\id\boxtimes\varphi_0}", from=1-5, to=1-7]
\arrow["{\rho_{F(X)}^\DD}"', from=1-5, to=2-5]
\arrow["{\varphi_{X,1}}", from=1-7, to=2-7]
\arrow["{F(\rho^\CC_X)}", from=2-7, to=2-5]
\end{tikzcd}\]

* (symmetry) For all $X,Y\in \CC$, the following diagram commutes.

\[\begin{tikzcd}
{F(X)\boxtimes F(Y)} && {F(Y)\boxtimes F(X)} \\
{F(X\otimes Y)} && {F(Y\otimes X)}
\arrow["{\sigma_{FX,FY}^\DD}", from=1-1, to=1-3]
\arrow["{\varphi_{X,Y}}"', from=1-1, to=2-1]
\arrow["{\varphi_{Y,X}}", from=1-3, to=2-3]
\arrow["{F(\sigma_{X,Y}^\CC)}"', from=2-1, to=2-3]
\end{tikzcd}\]

We write such a functor as a pair $(F,\varphi)$. When the comparison maps are isomorphisms, $(F,\varphi)$ is called a strong symmetric monoidal functor.

Given a pair of lax symmetric monoidal functors $(F,\varphi^F),(G,\varphi^G)\maps (\CC,\otimes, 1)\to (\DD,\boxtimes, e)$, a monoidal natural transformation from $(F,\varphi^F)$ to $(G,\varphi^G)$ is a natural transformation $\theta\maps F\Rightarrow G$ which also makes the following diagrams commute for all $X,Y\in \CC$.

\[\begin{tikzcd}
{F(X)\boxtimes F(Y)} && {G(X)\boxtimes G(Y)} && e && {F(1)} \\
{F(X\otimes Y)} && {G(X\otimes Y)} &&&& {G(1)}
\arrow["{\varphi^F_{X,Y}}"', from=1-1, to=2-1]
\arrow["{\varphi^G_{X,Y}}", from=1-3, to=2-3]
\arrow["{\theta_X\boxtimes\theta_Y}", from=1-1, to=1-3]
\arrow["{\theta_{X\otimes Y}}"', from=2-1, to=2-3]
\arrow["{\varphi^F_0}", from=1-5, to=1-7]
\arrow["{\theta_1}", from=1-7, to=2-7]
\arrow["{\varphi^G_0}"', from=1-5, to=2-7]
\end{tikzcd}\]


Counterfactual Reasoning for Bias Evaluation and Detection in a Fairness under Unawareness setting

A]Giandomenico Cornacchia0000-0001-5448-9970Corresponding Author. Email<EMAIL_ADDRESS>A]Vito Walter Anelli0000-0002-5567-4307
A]Fedelucio Narducci0000-0002-9255-3256
B]Azzurra Ragone0000-0002-3537-7663
A]Eugenio Di Sciascio0000-0002-5484-9945

[A]Polytechnic University of Bari
[B]Università degli Studi di Bari Aldo Moro

§ INTRODUCTION

§ RELATED WORK

§ PRELIMINARIES

§ METHODOLOGY

§ EXPERIMENTAL DESIGN

§ EXPERIMENTS AND DISCUSSION

§ CONCLUSION AND FUTURE WORK

This work was partially supported by the following projects: Secure Safe Apulia,
MISE CUP: I14E20000020001 CTEMT - Casa delle Tecnologie Emergenti Comune di Matera, CT_FINCONS_III, OVS Fashion Retail Reloaded, LUTECH DIGITALE 4.0, KOINÈ.

[1]
Sray Agarwal and Shashin Mishra, Responsible AI, Springer, 2021.

[2]
Ashwathy Ashokan and Christian Haas, `Fairness metrics and bias mitigation
strategies for rating predictions', Inf. Process. Manag., 58(5),
102646, (2021).

[3]
Mislav Balunovic, Anian Ruoss, and Martin T. Vechev, `Fair normalizing flows',
in The Tenth International Conference on Learning Representations,
ICLR 2022, Virtual Event, April 25-29, 2022. OpenReview.net, (2022).

[4]
P. J. Bickel, E. A. Hammel, and J. W. O'Connell, `Sex bias in graduate
admissions: Data from berkeley', Science, 187(4175),  398–404,

[5]
Léon Bottou, Jonas Peters, Joaquin Quiñonero-Candela, Denis X Charles,
D Max Chickering, Elon Portugaly, Dipankar Ray, Patrice Simard, and
Ed Snelson, `Counterfactual reasoning and learning systems: The example of
computational advertising.', Journal of Machine Learning Research, 14(11), (2013).

[6]
Consumer Financial Protection Bureau, `Using publicly available information to
proxy for unidentified race and ethnicity: a methodology and assessment',

[7]
Jiahao Chen, `Fair lending needs explainable models for responsible
recommendation', in FATREC'18 Proceedings of the Second Workshop on
Responsible Recommendation, Vancouver, British Columbia, Canada, (October 8

[8]
Jiahao Chen, Nathan Kallus, Xiaojie Mao, Geoffry Svacha, and Madeleine Udell,
`Fairness under unawareness: Assessing disparity when protected class is
unobserved', in FAT, pp. 339–348. ACM, (2019).

[9]
Sam Corbett-Davies, Emma Pierson, Avi Feller, Sharad Goel, and Aziz Huq,
`Algorithmic decision making and the cost of fairness', in KDD, pp.
797–806. ACM, (2017).

[10]
Giandomenico Cornacchia, Vito Walter Anelli, Giovanni Maria Biancofiore,
Fedelucio Narducci, Claudio Pomo, Azzurra Ragone, and Eugenio Di Sciascio,
`Auditing fairness under unawareness through counterfactual reasoning', Information Processing & Management, 60(2),  103224, (2023).

[11]
Giandomenico Cornacchia, Vito Walter Anelli, Fedelucio Narducci, Azzurra
Ragone, and Eugenio Di Sciascio, `Counterfactual reasoning for decision model
fairness assessment', in Companion Proceedings of the ACM Web Conference
2023, pp. 229–233, (2023).

[12]
Giandomenico Cornacchia, Vito Walter Anelli, Fedelucio Narducci, Azzurra
Ragone, and Eugenio Di Sciascio, `A general architecture for a trustworthy
creditworthiness-assessment platform in the financial domain', Annals of
Emerging Technologies in Computing (AETiC), 7(2), (2023).

[13]
Giandomenico Cornacchia, Fedelucio Narducci, and Azzurra Ragone, `A general
model for fair and explainable recommendation in the loan domain (short
paper)', in KaRS/ComplexRec@RecSys, volume 2960 of CEUR Workshop
Proceedings. CEUR-WS.org, (2021).

[14]
Kimberlé Williams Crenshaw, `Mapping the margins: Intersectionality,
identity politics, and violence against women of color', in The public
nature of private violence,  93–118, Routledge, (2013).

[15]
Sanjiv Das, Michele Donini, Jason Gelman, Kevin Haas, Mila Hardt, Jared
Katzman, Krishnaram Kenthapadi, Pedro Larroy, Pinar Yilmaz, and
Muhammad Bilal Zafar, `Fairness measures for machine learning in finance',
The Journal of Financial Data Science, 3(4),  33–64, (2021).

[16]
Michele Donini, Luca Oneto, Shai Ben-David, John Shawe-Taylor, and
Massimiliano Pontil, `Empirical risk minimization under fairness
constraints', in NeurIPS, pp. 2796–2806, (2018).

[17]
Julia Dressel and Hany Farid, `The accuracy, fairness, and limits of predicting
recidivism', Science Advances, 4(1),  eaao5580, (2018).

[18]
Cynthia Dwork, Moritz Hardt, Toniann Pitassi, Omer Reingold, and Richard S.
Zemel, `Fairness through awareness', in ITCS, pp. 214–226. ACM,

[19]
Marc Elliott, Allen Fremont, Peter Morrison, Philip Pantoja, and Nicole Lurie,
`A new method for estimating race/ethnicity and associated disparities where
administrative records lack self-reported race/ethnicity', Health
services research, 43,  1722–36, (05 2008).

[20]
Alessandro Fabris, Andrea Esuli, Alejandro Moreo, and Fabrizio Sebastiani,
`Measuring fairness under unawareness of sensitive attributes: A
quantification-based approach', J. Artif. Intell. Res., 76,
1117–1180, (2023).

[21]
Roberta Ferrario, `Counterfactual reasoning', in CONTEXT, volume 2116
of Lecture Notes in Computer Science, pp. 170–183. Springer, (2001).

[22]
Matthew L. Ginsberg, `Counterfactuals', Artif. Intell., 30(1),
35–79, (1986).

[23]
Moritz Hardt, Eric Price, and Nati Srebro, `Equality of opportunity in
supervised learning', in NIPS, pp. 3315–3323, (2016).

[24]
Hans Hofmann, UCI Statlog (German Credit Data) Data Set., UCI Machine
Learning Repository, 2000.

[25]
Max Hort and Federica Sarro, `Privileged and unprivileged groups: An empirical
study on the impact of the age attribute on fairness', in Proceedings of
the 2nd International Workshop on Equitable Data and Technology, FairWare
'22, p. 17–24, New York, NY, USA, (2022). Association for Computing

[26]
Amir-Hossein Karimi, Gilles Barthe, Borja Balle, and Isabel Valera,
`Model-agnostic counterfactual explanations for consequential decisions', in
AISTATS, volume 108 of Proceedings of Machine Learning
Research, pp. 895–905. PMLR, (2020).

[27]
Ronny Kohavi and Barry Becker, UCI Adult Data Set., UCI Machine Learning
Repository, May 1996.

[28]
Matt J. Kusner, Joshua R. Loftus, Chris Russell, and Ricardo Silva,
`Counterfactual fairness', in NIPS, pp. 4066–4076, (2017).

[29]
Ninareh Mehrabi, Fred Morstatter, Nripsuta Saxena, Kristina Lerman, and Aram
Galstyan, `A survey on bias and fairness in machine learning', ACM
Computing Surveys (CSUR), 54(6),  1–35, (2021).

[30]
Tim Miller, `Explanation in artificial intelligence: Insights from the social
sciences', Artif. Intell., 267,  1–38, (2019).

[31]
Alan Mishler, Edward H. Kennedy, and Alexandra Chouldechova, `Fairness in risk
assessment instruments: Post-processing to achieve counterfactual equalized
odds', in FAccT, pp. 386–400. ACM, (2021).

[32]
Ramaravind K Mothilal, Amit Sharma, and Chenhao Tan, `Explaining machine
learning classifiers through diverse counterfactual explanations', in Proceedings of the 2020 Conference on Fairness, Accountability, and
Transparency, pp. 607–617, (2020).

[33]
Ramaravind K. Mothilal, Amit Sharma, and Chenhao Tan, `Explaining machine
learning classifiers through diverse counterfactual explanations', in Proceedings of the 2020 Conference on Fairness, Accountability, and
Transparency, FAT* '20, p. 607–617, New York, NY, USA, (2020). Association
for Computing Machinery.

[34]
Judea Pearl, `Causation, action and counterfactuals', in ECAI, pp.
826–828. John Wiley and Sons, Chichester, (1994).

[35]
Dino Pedreschi, Salvatore Ruggieri, and Franco Turini, `Discrimination-aware
data mining', in KDD, pp. 560–568. ACM, (2008).

[36]
Michael Redmond, UCI Communities and Crime Data Set, UCI Machine Learning
Repository, 1995.

[37]
Boris Ruf and Marcin Detyniecki, `Active fairness instead of unawareness', CoRR, abs/2009.06251, (2020).

[38]
Yongkai Wu, Lu Zhang, and Xintao Wu, `Counterfactual fairness:
Unidentification, bound and algorithm', in IJCAI, pp. 1438–1444.
ijcai.org, (2019).

[39]
Samuel Yeom, Anupam Datta, and Matt Fredrikson, `Hunting for discriminatory
proxies in linear regression models', Advances in Neural Information
Processing Systems, 31, (2018).

[40]
Muhammad Bilal Zafar, Isabel Valera, Manuel Gomez Rogriguez, and Krishna P.
Gummadi, `Fairness Constraints: Mechanisms for Fair Classification', in
Proceedings of the 20th International Conference on Artificial
Intelligence and Statistics, eds., Aarti Singh and Jerry Zhu, volume 54 of
Proceedings of Machine Learning Research, pp. 962–970. PMLR, (20–22
Apr 2017).

[41]
Brian Hu Zhang, Blake Lemoine, and Margaret Mitchell, `Mitigating unwanted
biases with adversarial learning', in Proceedings of the 2018 AAAI/ACM
Conference on AI, Ethics, and Society, AIES '18, p. 335–340, New York, NY,
USA, (2018). Association for Computing Machinery.

[42]
Ziwei Zhu, Xia Hu, and James Caverlee, `Fairness-aware tensor-based
recommendation', in CIKM, pp. 1153–1162. ACM, (2018).

§ APPENDIX