#  AntibodyFlow: Normalizing Flow Model for Designing Antibody
Complementarity-Determining Regions

Bohao Xu1 Yanbo Wang2 Wenyu Chen2 and Shimin Shan3

(1 School of Electrical Engineering and Computer Science, Nanjing University
of Aeronautics and Astronautics, NanJing
2 School of Computer Science and Technology, North University of China,
Taiyuan
3 School of Semiconductor and Physics, North University of China, Taiyuan
)

###### Abstract

Therapeutic antibodies have been extensively studied in drug discovery and
development in the past decades. Antibodies are specialized protective
proteins that bind to antigens in a lock-to-key manner. The binding
strength/affinity between an antibody and a specific antigen is heavily
determined by the complementarity-determining regions (CDRs) on the
antibodies. Existing machine learning methods cast in silico development of
CDRs as either sequence or 3D graph (with a single chain) generation tasks and
have achieved initial success. However, with CDR loops having specific
geometry shapes, learning the 3D geometric structures of CDRs remains a
challenge. To address this issue, we propose AntibodyFlow, a 3D flow model to
design antibody CDR loops. Specifically, AntibodyFlow first constructs the
distance matrix, then predicts amino acids conditioned on the distance matrix.
Also, AntibodyFlow conducts constraint learning and constrained generation to
ensure valid 3D structures. Experimental results indicate that AntibodyFlow
outperforms the best baseline consistently with up to 16.0% relative
improvement in validity rate and 24.3% relative reduction in geometric graph
level error (root mean square deviation, RMSD).

##  1 Introduction

Therapeutic antibodies have been extensively studied in drug discovery and
development in the past decades [1]. Compared to small-molecule drugs,
therapeutic antibodies have many advantages such as minimal adverse effects,
and can bind to many previously “undruggable” targets due to different
biochemical mechanisms [2]. The majority of the affinity and specificity of
antibodies is modulated by a set of binding loops called the Complementarity
Determining Regions (CDRs) found on the variable domain of each of the two
chains. One critical task for antibody discovery is to design the appropriate
CDR [3, 4]. Traditionally, antibody sequences are discovered either in vivo
using animal immunization or by in vitro affinity selection of candidates from
large synthetic libraries of antibody sequences [5]. However, these approaches
are time- and labor-intensive and unable to efficiently explore the large
antibody space since the design process heavily relies on existing human
knowledge [3].

In silico antibody CDR design methods were developed to tackle this challenge.
Earlier methods are primarily based on heuristic searching and require
extensive wet-lab validation [6, 7, 8]. In recent years, several sequence/3D
graph generation methods have been proposed for antibody CDR loop discovery,
e.g., [9, 10, 11]. However, several challenges remain:

  1. 1.

Most of the existing methods generate amino acids sequentially, which is hard
to model long-term dependency and impose global constraints on the geometric
shape of CDR loops.

  2. 2.

Rotating and translating the 3D graph would not change the underlying 3D graph
and it is challenging to directly generate rotation- and translation-invariant
3D graphs [12].

To address these challenges, we propose AntibodyFlow, a 3D flow model that
generates antibody CDR loops in two phases: (i) build the geometric structure
of the CDR loop with spatial constraints; and (ii) predict the amino acids as
the graph nodes. The main contributions are summarized as follows.

  1. 1.

AntibodyFlow generates amino acid sequence and pairwise distance matrix
(between amino acids) in a one-hot manner. Unlike the sequential manner, the
flow model with one-hot generation is more flexible to impute the global
constraints of geometric shape. Concretely, AntibodyFlow designs constraint
learning to learn the constraints and constrains 3D coordinates generation
leveraging the constraints.

  2. 2.

To circumvent the direct generation of 3D graphs, we build a pairwise distance
matrix that measures Euclidean distances between amino acids. The distance
matrix is invariant w.r.t. rotation and translation on a 3D graph. After that,
we generate 3D coordinates based on both distance matrix and geometric
constraints.

  3. 3.

Experiments confirmed AntibodyFlow outperforms the best baseline consistently
with up to 16.0% relative improvement in validity rate and 24.3% relative
reduction in geometric graph level error (root mean square deviation, RMSD).

##  2 Method

###  2.1 Problem Formulation

This section describes the antibody CDR loop and formulates the antibody CDR
design problem. For ease of exposition, we list the important mathematical
notations in Table 3 in Appendix.

Amino acid. A protein (including antibodies) is an amino acid sequence folded
into a three-dimensional structure. Specifically, amino acids are small
organic molecules that consist of an $\alpha$ (central) carbon atom linked to
an amino group, a carboxyl group, a hydrogen atom, and a variable component
called a side chain. Multiple amino acids are linked together within a protein
by peptide bonds, thereby forming a long chain. Amino acids are the basic
building blocks of a protein/antibody. The set of all the natural amino acids
(denoted $V$) consists of 20 natural amino acids (i.e., $|V|=20$), e.g.,
Alanine, Cysteine, Glutamine, etc.

Figure 1: Data representation. An antibody is a special kind of protein with a
symmetric Y shape, each half of the symmetric unit has two chains: a heavy
chain (H) and a light chain (L). In total, there are four chains, two
identical H/L chains. The majority of the binding affinity (to specific
antigens) is modulated by a set of binding loops called the Complementarity
Determining Regions (CDRs) found on the variable domain of each of the H and L
chains. There are 6 CDR loops on each half of the antibody, L1, L2, L3 on the
light chain, and H1, H2, H3 on the heavy chain. We show the H3 loop of the
antibody whose protein data bank (PDB) ID is “5iwl” as an example. The CDR
loop in the 3D geometry graph contains 7 amino acids (GGYRAMD) and their
coordinates. The sequence is represented as a binary matrix
$\mathbf{S}\in\\{0,1\\}^{7\times 20}$ (20 natural amino acids, each row is a
one-hot vector), and the geometry information is represented as a pairwise
distance matrix $\mathbf{D}\in\mathbb{R}^{7\times 7}_{+}$.

CDR loops. An antibody is a special kind of protein with a symmetric Y shape;
each half of the symmetric unit has two chains: a heavy chain (H) and a light
(L) chain, as shown in Figure 1. There are six CDR loops, L1, L2, L3 on the
light chain and H1, H2, H3 on the heavy chain [3, 4]. CDR is complicated to
model computationally, causing a significant hurdle for in silico development
of antibody biotherapeutics. Compared with L1, L2, and L3 loops, the H1, H2,
and H3 loop in the CDR of an antibody plays a critical role in its binding
ability to potential antigens [1, 13]. Thus, this paper restricts our
attention to modeling H1, H2, and H3 loops separately.

Data representation. A 3D CDR loop (H1, or H2, or H3) can be represented as
$\mathbf{S}$ amino acid sequence and their coordinates $\mathbf{G})$. Given a
CDR loop of $N$ amino acids, $\mathbf{S}\in\\{0,1\\}^{N\times|V|}$ is a binary
matrix, the $i$-th row of $\mathbf{S}$ is a one-hot vector that indicates the
category of the $i$-th amino acid;
$\mathbf{G}=[\mathbf{g}_{1},\cdots,\mathbf{g}_{N}]^{\top}\in\mathbb{R}^{N\times
3}$ represents the coordinates for all the amino acids (coordinates of all
$\alpha$ atoms). Since rotating or translating $\mathbf{G}$ would not change
the underlying 3D graph, it is challenging to generate such a 3D graph
$\mathbf{G}$ [14]. To circumvent the direct generation of 3D amino acids’
coordinates $\mathbf{G}$, we instead leverage the pairwise distance matrix
$\mathbf{D}\in\mathbb{R}_{+}^{N\times N}$, which describes the Euclidean
distance between the amino acids. That is, the $i,j$-th element of
$\mathbf{D}$ satisfies
$\mathbf{D}_{i,j}=||\mathbf{g}_{i}-\mathbf{g}_{j}||_{2}$. The CDR loop can be
represented as $\mathbf{x}=(\mathbf{D},\mathbf{S})$. An illustration of data
representation is given in Figure 1.

Now, we show that distance matrix $\mathbf{D}$ is invariant with respect to
translation and rotation on the 3D graph whose coordinates are $\mathbf{G}$.
That is, rotating and translating $\mathbf{G}$ would not change $\mathbf{D}$.

######  Corollary 1.

Suppose $\mathbf{G}^{\prime}$ are the 3D coordinates of the graph by rotating
and translating the original graph whose coordinates are $\mathbf{G}$,
$\mathbf{D}$ and $\mathbf{D}^{\prime}$ are the distance matrix of $\mathbf{G}$
and $\mathbf{G}^{\prime}$, respectively, then we have
$\mathbf{D}=\mathbf{D}^{\prime}$. The proof is given in the Appendix.

Validity constraints. We define the validity of the generated loop based on
empirical domain knowledge about CDR loops [15, 16, 17]. Specifically, we
define the validity of a generated 3D CDR loop when it satisfies the following
two constraints:

  1. 1.

Fixed bond length. The amino acid bonds have relatively fixed length, i.e.,
the distance between connected amino acids are relatively fixed [15]:

$\eta_{1}\leq\mathbf{D}_{i,i+1}\leq\eta_{2},\ \text{for}\ i=1,\cdots,N-1.$ (1)

  2. 2.

Open loop. The shape of CDR is an open loop, where the distance between the
first and the last amino acids is within a specific range [16, 17],

$\epsilon_{1}\leq\mathbf{D}_{1,N}\leq\epsilon_{2}.$ (2)

The setup of $\eta_{1},\eta_{2},\epsilon_{1},\epsilon_{2}$ is based on domain
knowledge [16, 17] and empirical validation. For example, for CDR H3 loop, we
set $\eta_{1}=3.71\mathring{A}$, $\eta_{2}=3.88\mathring{A}$,
$\epsilon_{1}=6.5\mathring{A},\epsilon_{2}=8.5\mathring{A}$, where Angstrom
$\mathring{A}$ is a metric unit of length equal to $10^{-10}$ m. Setups for H1
and H2 loops are given in Section H in Appendix.

###  2.2 Background and Overview

Flow model background. The objective of flow model is to learn an invertible
model $f_{\theta}$ (with learnable parameter $\theta$) that maps the input
data ${\mathbf{x}}\in\mathbb{R}^{d}$ into the latent variable
${\mathbf{z}}\in\mathbb{R}^{d}$, i.e.,
${\mathbf{z}}=f_{\theta}({\mathbf{x}})$. $f_{\theta}$ is invertible, i.e.,
${\mathbf{x}}=f^{-1}_{\theta}({\mathbf{z}})$. The latent variable has the same
size as the input data. We suppose the latent variable ${\mathbf{z}}$ is drawn
from a known prior distribution (e.g., isotropic Gaussian distribution), i.e.,
${\mathbf{z}}\sim p_{{\mathbf{z}}}({\mathbf{z}})$. Using the change-of-
variable formula, the log-probability of ${\mathbf{x}}$ becomes

$\log p_{{\mathbf{x}}}({\mathbf{x}})=\log
p_{\mathbf{z}}({\mathbf{z}})+\log\bigg{(}\bigg{|}\det\frac{\partial{\mathbf{z}}}{\partial{\mathbf{x}}}\bigg{|}\bigg{)},$
(3)

where $\frac{\partial{\mathbf{z}}}{\partial{\mathbf{x}}}$ is the Jacobian of
$f_{\theta}$ at $\mathbf{x}$, $\det$ is the determinant of a square matrix.
Recall input data $\mathbf{x}$ is
$\mathbf{x}=(\mathbf{D},\mathbf{S}),\mathbf{D}\in\mathbb{R}^{N\times
N},\mathbf{S}\in\\{0,1\\}^{N\times|V|}$. The design of $f$ follows two
strategies to learn complex transformations: (i) “easy logarithm of
determinant of the Jacobian”
($\log\big{|}\det\frac{\partial\mathbf{z}}{\partial\mathbf{x}}\big{|}$) and
(ii) “easy inverse” ($f^{-1}$) [18].

Background of affine coupling. Affine coupling layer is a commonly-used neural
architecture in flow model to learn complex transformation, which naturally
meets the two criteria above [18, 19, 20]. Specifically, in affine coupling
layer, the input $\mathbf{I}\in\mathbb{R}^{P\times Q}$ is decomposed into
$\mathbf{I}=(\mathbf{I}_{1},\mathbf{I}_{2}),\mathbf{I}_{1}\in\mathbb{R}^{P_{1}\times
Q},\mathbf{I}_{2}\in\mathbb{R}^{(P-P_{1})\times Q}$, where
$P_{1}=\lceil{\frac{P}{2}}\rceil$ (ceiling function) or
$P-\lceil{\frac{P}{2}}\rceil$. The mapping is

$\mathbf{O}_{1}=\mathbf{I}_{1},\ \
\mathbf{O}_{2}=\mathbf{I}_{2}\odot\text{sigmoid}(R(\mathbf{I}_{1}))+T(\mathbf{I}_{1}),$
(4)

where $\odot$ denotes element-wise multiplication. Outputs
$\mathbf{O}_{1}\in\mathbb{R}^{P_{1}\times
Q},\mathbf{O}_{2}\in\mathbb{R}^{(P-P_{1})\times Q}$ has same shape with the
inputs $\mathbf{I}_{1}\in\mathbb{R}^{P_{1}\times
Q},\mathbf{I}_{2}\in\mathbb{R}^{(P-P_{1})\times Q}$, respectively.
$R(\cdot),T(\cdot):\mathbb{R}^{P_{1}\times
Q}\xrightarrow[]{}\mathbb{R}^{(P-P_{1})\times Q}$ share the same neural
architectures but different parameters. The sigmoid function is applied
element-wisely. The logarithm of the Jacobian determinant is easy to evaluate,

$\displaystyle\log\bigg{|}\det\frac{\partial\mathbf{O}}{\partial\mathbf{I}}\bigg{|}=\log\bigg{|}\det\begin{pmatrix}\mathbb{I}&0\\\
0&\text{sigmoid}(R(\mathbf{I}_{1}))\\\
\end{pmatrix}\bigg{|}=\sum\nolimits_{j=1}^{(P-P_{1})Q}\log\text{sigmoid}(R(\mathbf{I}_{1}))_{j}.$
(5)

where $\mathbb{I}$ is the identity matrix. The inverse mapping is

$\displaystyle\mathbf{I}_{1}=\mathbf{O}_{1},\ \
\mathbf{I}_{2}=\big{(}\mathbf{O}_{2}-T(\mathbf{I}_{1})\big{)}/\text{sigmoid}\big{(}R(\mathbf{I}_{1})\big{)}.$
(6)

Multiple affine coupling layers are usually concatenated (e.g.,
$\mathcal{F}_{1}$, $\cdots$, $\mathcal{F}_{K}$ in Eq.7). Consecutive affine
coupling layers swap the elements in $\mathbf{I}_{1}$ and $\mathbf{I}_{2}$ so
that every element occurs in $\mathbf{I}_{1}$ and $\mathbf{I}_{2}$ in turn.

Overview of AntibodyFlow. AntibodyFlow generates 3D CDR loop $\mathbf{x}$ in
two phases. The first phase generates the distance matrix $\mathbf{D}$. In the
second phase, based on $\mathbf{D}$, we predict the amino acids for all the
nodes, denoted $\mathbf{S}$. Specifically, the latent variable ${\mathbf{z}}$
is decomposed as
${\mathbf{z}}=({\mathbf{z}}_{\mathbf{D}},{\mathbf{z}}_{\mathbf{S}|\mathbf{D}})$,
the invertible model $f=(f_{\mathbf{D}},f_{\mathbf{S}|\mathbf{D}})$111$f$
should be formally written as $f_{\theta}$. We omit $\theta$ for simplicity.
has two components. The two latent variables are
${\mathbf{z}}_{\mathbf{D}}=f_{\mathbf{D}}(\mathbf{D}),\mathbf{z}_{\mathbf{D}}\in\mathbb{R}^{N\times
N}$, and
${\mathbf{z}}_{\mathbf{S}|\mathbf{D}}=f_{\mathbf{S}|\mathbf{D}}(\mathbf{S}|\mathbf{D}),\mathbf{z}_{\mathbf{S}|\mathbf{D}}\in\mathbb{R}^{N\times|V|}$.
The Jacobian of ${\mathbf{z}}=f({\mathbf{x}})$ is decomposed as

$\frac{\partial{\mathbf{z}}}{\partial{\mathbf{x}}}=\frac{\partial({\mathbf{z}}_{\mathbf{D}},{\mathbf{z}}_{\mathbf{S}|\mathbf{D}})}{\partial(\mathbf{D},\mathbf{S})}=\begin{bmatrix}\frac{\partial{\mathbf{z}}_{\mathbf{D}}}{\partial\mathbf{D}}&\frac{\partial{\mathbf{z}}_{\mathbf{D}}}{\partial\mathbf{S}}\\\
\frac{\partial{\mathbf{z}}_{\mathbf{S}|\mathbf{D}}}{\partial\mathbf{D}}&\frac{\partial{\mathbf{z}}_{\mathbf{S}|\mathbf{D}}}{\partial\mathbf{S}}\\\
\end{bmatrix}=\begin{bmatrix}\frac{\partial{\mathbf{z}}_{\mathbf{D}}}{\partial\mathbf{D}}&0\\\
\frac{\partial{\mathbf{z}}_{\mathbf{S}|\mathbf{D}}}{\partial\mathbf{D}}&\frac{\partial{\mathbf{z}}_{\mathbf{S}|\mathbf{D}}}{\partial\mathbf{S}}\\\
\end{bmatrix}$

The logarithm Jacobian determinant of $f_{\theta}$ can be decomposed as
$\log\big{|}\det\frac{\partial\mathbf{z}}{\partial\mathbf{x}}\big{|}=\log\big{|}\det\frac{\partial({\mathbf{z}}_{\mathbf{D}},{\mathbf{z}}_{\mathbf{S}|\mathbf{D}})}{\partial(\mathbf{D},\mathbf{S})}\big{|}=\log\big{|}\det\frac{\partial{\mathbf{z}}_{\mathbf{D}}}{\partial\mathbf{D}}\big{|}+\log\big{|}\det\frac{\partial{\mathbf{z}}_{\mathbf{S}|\mathbf{D}}}{\partial\mathbf{S}}\big{|}$.
The flow model is split into two parts: (1) distance matrix flow (Section
2.3); (2) conditional amino acid flow (conditioned on distance matrix)
(Section 2.4). Section 2.5 discusses how flow model learns constraints;
Section 2.6 describes how to generate 3D coordinates based on both distance
matrix and constraints. Section D (Appendix) describes the learning and
inference process. The pipeline is illustrated in Fig. 3. We slightly abuse
the mathematical notations for simplicity.

###  2.3 Distance Matrix Flow

This section elaborates on distance matrix flow $f_{\mathbf{D}}(\mathbf{D})$.
We illustrate the neural architecture in Figure 3. To build invertible
mappings, we concatenate $K$ affine coupling layers (Eq. 4, denoted
$\mathcal{F}_{1},\cdots,\mathcal{F}_{K}$), which is a popular technique in
flow model [19, 20],

$f_{\mathbf{D}}=\mathcal{F}_{1}\circ\mathcal{F}_{2}\circ\cdots\circ\mathcal{F}_{K}.$
(7)

The logarithm of Jacobian determinant of $f_{\mathbf{D}}$ is

$\log\bigg{|}\det\frac{\partial
f_{\mathbf{D}}}{\partial\mathbf{D}}\bigg{|}=\sum\nolimits_{k=1}^{K}\log\bigg{|}\det\frac{\partial\mathcal{F}_{k}}{\partial
Y_{k}}\bigg{|}.$ (8)

where $Y_{k}$ is the input of $k$-th affine coupling layer, and also the
output of ($k$-1)-th affine coupling layer,
$\big{|}\det\frac{\partial\mathcal{F}_{k}}{\partial Y_{k}}\big{|}$ is given in
Eq. (5). The inverse flow $f_{\mathbf{D}}^{-1}$ reverses the order to apply
the $K$ affine coupling layers (Eq. 7):

$f^{-1}_{\mathbf{D}}=\mathcal{F}_{K}^{-1}\circ\mathcal{F}_{K-1}^{-1}\circ\cdots\circ\mathcal{F}_{1}^{-1}.$
(9)

The inverse of affine coupling ($\mathcal{F}^{-1}_{k}$) is given in Eq. (6).
The mapping of an affine coupling layer is given in Eq. 4. As mentioned,
scale/transformation function $R(\cdot)/T(\cdot)$ (Eq. 4) have the same neural
architecture: to capture the local relationship between the distance matrix,
we leverage a two-dimensional convolutional neural network (CNN2D) followed by
a batch normalization (BN) layer and ReLU activation, shown in Fig 3.

###  2.4 Conditional Amino Acid Flow

This section describes the conditional amino acid flow
$f_{\mathbf{S}|\mathbf{D}}(\mathbf{S}|\mathbf{D})$, which shares similarity
with the distance matrix flow in Section 2.3. The main difference is the
neural architecture of $R(\cdot)/T(\cdot)$, as shown in Figure 3.

Conditional amino acid flow model $f_{\mathbf{S}|\mathbf{D}}$ concatenates $L$
affine coupling layers, denoted $\mathcal{H}_{1},\cdots,\mathcal{H}_{L}$

$f_{\mathbf{S}|\mathbf{D}}=\mathcal{H}_{1}\circ\mathcal{H}_{2}\circ\cdots\circ\mathcal{H}_{L}.$
(10)

Similar to Eq. (8) and (9), the logarithm of Jacobian determinant of
$f_{\mathbf{S}|\mathbf{D}}$ and the inverse are

$\displaystyle\log\bigg{|}\det\frac{\partial
f_{\mathbf{S}|\mathbf{D}}}{\partial\mathbf{S}}\bigg{|}=\sum\nolimits_{l=1}^{L}\log\bigg{|}\det\frac{\partial\mathcal{H}_{l}}{\partial
Y_{l}}\bigg{|},\ \ \ \
f^{-1}_{\mathbf{S}|\mathbf{D}}=\mathcal{H}_{L}^{-1}\circ\mathcal{H}_{L-1}^{-1}\circ\cdots\circ\mathcal{H}_{1}^{-1}.$
(11)

where $Y_{l}$ is the input of $l$-th affine coupling layer and also the output
of ($l$-1)-th layer. For $\mathcal{H}_{l}$ ($l$=1,$\cdots$,$L$), the mapping
of an affine coupling layer is
$\mathbf{O}_{1}=\mathbf{I}_{1},\mathbf{O}_{2}=\mathbf{I}_{2}\odot\text{sigmoid}\big{(}R(\mathbf{I}_{1}|\mathbf{D})\big{)}+T(\mathbf{I}_{1}|\mathbf{D})$
(Eq. 4). $\big{|}\det\frac{\partial\mathcal{H}_{l}}{\partial Y_{l}}\big{|}$
and $\mathcal{H}_{1}^{-1}$ are in Eq. (5) and (6), respectively. Consecutive
affine coupling layers ($\mathcal{H}_{1},\cdots,\mathcal{H}_{K}$) swap the
elements in $\mathbf{I}_{1}$ and $\mathbf{I}_{2}$.

Scale/transformation functions $R(\cdot)/T(\cdot)$ (Eq. 4) have same
architecture. The architecture consists of (1) wGNN (weighted-distance graph
neural network, Eq.12); (2) Batch Normalization (BN), and (3) ReLU activation.
The traditional graph neural networks (GNNs), e.g., graph convolutional
network [21], graph attention network [22], pass the message between nodes
through binary adjacency matrix. In contrast, in our scenario, the CDR loop is
a single chain and is folded into the 3D structure. To characterize the
interaction between nodes, we design weighted-distance graph neural network to
incorporate weighted edge, denoted wGNN,

$\text{wGNN}(\mathbf{I}_{1}|\mathbf{D})=\tilde{\mathbf{D}}(\mathbf{M}(\mathbf{I}_{1})\odot\mathbf{I})\mathbf{W},$
(12)

where $\mathbf{W}\in\mathbb{R}^{|V|\times|V|}$ are the learnable weight
matrix; $\mathbf{M}(\mathbf{I}_{1})\in\\{0,1\\}^{N\times|V|}$ is a binary mask
to select a partition $\mathbf{I}_{1}$ from
$\mathbf{I}\in\mathbb{R}^{N\times|V|}$;
$\tilde{\mathbf{D}}\in\mathbb{R}^{N\times N}$ is weighted edge matrix; the
$(i,j)$-th element $\tilde{\mathbf{D}}_{i,j}=\exp(-\xi\mathbf{D}_{i,j})$
measures the interaction strength between $i$-th and $j$-th nodes, $\xi>0$ is
the hyperparameter. Intuitively, when the $i$-th and $j$-th nodes are closer,
the interaction strength is stronger. In sum, $\tilde{\mathbf{D}}$ replaces
the binary adjacency matrix in the vanilla graph neural network.

###  2.5 Constraint Learning

Constraint learning encourages the flow model to learn the validity
constraints of the CDR loop. To proceed, we design a constraint loss function,
a function of the generated distance matrix $\mathbf{D}$. It penalizes the
distance matrix that does not satisfy the validity constraints. However, the
constraints in Eq. (1) and (2) are hard constraints and not differentiable
w.r.t. the distance matrix. To bridge the gap, we design the following
constraint surrogate function $H(y;a,b,\delta)$ that transforms the hard
constraints into soft constraints.

The constraint surrogate function $H(y;a,b,\delta)$ is a function of $y$ (a
scalar), where $a,b,\delta$ are all scalars and satisfy $a<b$ and $\delta>0$.
$H$ is designed to encourage the variable $y$ to be between scalar $a$ and
$b$. To make $H$ differentiable everywhere and enable the continuous
optimization process, we borrow the idea from Huber loss [23]. The basic idea
of Huber loss comes from the fact that square loss is sensitive to outliers
but differentiable everywhere, whereas L1 loss is more robust to outliers but
not differentiable at 0. Combining the advantage of both, the Huber loss
function is quadratic for small values (square loss), and linear for large
values (L1 loss). In our scenario, as long as the scalar $y$ is between $a$
and $b$, the constraints are satisfied with the loss equal to 0. Thus, we
design $H(y;a,b,\delta)$ as

$\displaystyle H(y;a,b,\delta)=\left\\{\begin{aligned}
&\delta(-y+a-\frac{1}{2}\delta),&&y\leq a-\delta\\\
&\frac{1}{2}(y-a)^{2},&&a-\delta<y\leq a,\\\ &0,&&a<y\leq b,\\\
&\frac{1}{2}(y-b)^{2},&&b<y\leq b+\delta,\\\
&\delta(y-b-\frac{1}{2}\delta),&&b+\delta<y.\\\ \end{aligned}\right.$ (13)

It has equal values and slopes of the different sections at the four
connection points ($a-\delta,a,b,b+\delta$), so it is continuous and
differentiable everywhere. When $a=b$, $H$ loss function reduces to Huber loss
function. Thus, the hard constraints in Eq. (1) and (2) can be transformed
into surrogate functions as follows

$\displaystyle
h_{1}(\mathbf{D})=\sum\nolimits_{i=1}^{N-1}H(\mathbf{D}_{i,i+1};\eta_{1},\eta_{2},\eta_{3}),\
\ \ \ \
h_{2}(\mathbf{D})=H(\mathbf{D}_{1,N};\epsilon_{1},\epsilon_{2},\epsilon_{3}).$
(14)

where $\eta_{1}/\eta_{2}$ and $\epsilon_{1}/\epsilon_{2}$ are defined in Eq.
(1) and (2), respectively. We also have $\eta_{3}=\eta_{2}-\eta_{1}$ and
$\epsilon_{3}=\epsilon_{2}-\epsilon_{1}$.

Algorithm 1 AntibodyFlow

1: Input: data ${\mathbf{x}}=\\{\mathbf{S},\mathbf{D}\\}$; Training iteration
$T$.

2: Output: $\widetilde{{X}^{(1)}},\widetilde{{X}^{(2)}},\cdots$.

3: #### Learning ####

4: for $t=1,2,\cdots,T$ do

5: Maximize log-likelihood $\log p_{\mathbf{x}}(\mathbf{x})$(Eq.3)(Sec2.3-2.4)

6: Minimize constraint loss $\widetilde{\mathcal{L}_{\text{con}}}$ (Eq. 19)
(Sec 2.5).

7: end for

8: #### Inference ####

9: for $i=1,2,\cdots,$ do

10: Draw latent variable ${\mathbf{z}}\sim p_{\mathbf{z}}({\mathbf{z}})$, and
${\mathbf{z}}=({\mathbf{z}}_{\mathbf{D}},{\mathbf{z}}_{\mathbf{S}|\mathbf{D}})$.

11: Generate distance matrix:
$\mathbf{D}=f_{\mathbf{D}}^{-1}({\mathbf{z}}_{\mathbf{D}})$ (Eq. 9).

12: Generate amino acid matrix:
$\mathbf{S}=f_{\mathbf{S}|\mathbf{D}}^{-1}({\mathbf{z}}_{\mathbf{S}|\mathbf{D}};\mathbf{D})$
(Eq. 11).

13: Constrained 3D coordinates generation: $\mathbf{G}$ (Eq. 20) (Sec 2.6).

14: $\widetilde{{X}^{(i)}}=(\mathbf{S},\mathbf{G})$.

15: end for

$h_{1}(\cdot)$ and $h_{2}(\cdot)$ consider $\mathbf{D}_{1,N}$ and
$\mathbf{D}_{i,i+1}$ ($i=1,\cdots,N-1$), respectively. On the other hand,
based on the validity constraint (Eq.1) and triangle inequality, we induce
some constraints on $\\{\mathbf{D}_{i,j}\\}$ ($i,j\in\\{1,\cdots,N\\}$) as
follows,

$\displaystyle|\mathbf{D}_{i,j}-\mathbf{D}_{(i+1),j}|\leq\mathbf{D}_{i,(i+1)}\leq\eta_{2},\
\ \
|\mathbf{D}_{i,j}-\mathbf{D}_{i,(j+1)}|\leq\mathbf{D}_{j,(j+1)}\leq\eta_{2},$
(15) $\displaystyle\text{for}\ i,j\in\\{1,2,\cdots,N-1\\},$

where $\eta_{2}$ is defined in Eq. (1). It indicates that the elements of
$\mathbf{D}$ are smooth, we design the following constraint loss to encourage
the smoothness,

$\displaystyle
h_{3}(\mathbf{D})=\sum_{i=1}^{N-1}\sum_{j=1}^{N-1}(\mathbf{D}_{i,j}-\mathbf{D}_{(i+1),j})^{2}+(\mathbf{D}_{i,j}-\mathbf{D}_{i,(j+1)})^{2}.$
(16)

We design constraint loss function as follows, which is a function of distance
matrix $\mathbf{D}$,

$\displaystyle\
h(\mathbf{D})=\xi_{1}h_{1}(\mathbf{D})+\xi_{2}h_{2}(\mathbf{D})+\xi_{3}h_{3}(\mathbf{D}),$
(17)

where $\xi_{1},\xi_{2},\xi_{3}>0$ are the hyperparameters that control the
weights of various terms. When generating new distance matrix, we draw latent
variable ${\mathbf{z}}_{\mathbf{D}}$ from isotropic Gaussian distribution
$p_{{\mathbf{z}}_{\mathbf{D}}}(\cdot)$, the learning objective is designed as
expectation of $h(\mathbf{D})$ under the distribution of
$\mathbf{z}_{\mathbf{D}}$

$\displaystyle\mathcal{L}_{\text{con}}(\theta)=\mathbb{E}_{\big{[}{\mathbf{z}}_{\mathbf{D}}\sim
p_{{\mathbf{z}}_{\mathbf{D}}}(\cdot),\mathbf{D}=(f_{\theta})_{\mathbf{D}}^{-1}({\mathbf{z}}_{\mathbf{D}})\big{]}}h(\mathbf{D}).$
(18)

We use Monte Carlo sampling to get a noisy but unbiased estimate of the
learning objective,

$\displaystyle\widetilde{\mathcal{L}_{\text{con}}(\theta)}=\frac{1}{M}\sum\nolimits_{i=1}^{M}h(\mathbf{D}^{(i)}),\
\text{where}\
{\mathbf{z}}^{(i)}_{\mathbf{D}}\overset{\text{i.i.d.}}{\sim}p_{{\mathbf{z}}_{\mathbf{D}}}(\cdot),$
(19) $\displaystyle i=1,\cdots,M,\ \ \ \ \ \
\mathbf{D}^{(i)}=(f_{\theta})_{\mathbf{D}}^{-1}({\mathbf{z}}^{(i)}_{\mathbf{D}}),$

where $M$ is the number of Monte Carlo samples. Gradient descent is leveraged
to minimize the noisy objective
$\widetilde{\mathcal{L}_{\text{con}}(\theta)}$. We alternatively optimize log-
likelihood in Eq. (3) and constraint loss in Eq. (17). Constraint learning is
only conducted during learning and not utilized in inference.

###  2.6 Constrained 3D Coordinates Generation

This section discusses how to transform the generated distance matrix
$\mathbf{D}$ into 3D coordinates $\mathbf{G}$ that satisfy validity
constraints simultaneously. Specifically, given the generated distance matrix
$\mathbf{D}\in\mathbb{R}^{N\times N}$, we build the 3D graph via optimizing
the coordinates of all the amino acids, i.e.,
$\mathbf{G}=[\mathbf{g}_{1},\mathbf{g}_{2},\cdots,\mathbf{g}_{N}]^{\top}\in\mathbb{R}^{N\times
3}$. We optimize the amino acids’ coordinates by minimizing the following
objective under the validity constraints given in Eq. (1) and (2):

$\displaystyle\underset{\mathbf{g}_{1},\cdots,\mathbf{g}_{N}\in\mathbb{R}^{3}}{\arg\min}\sum\nolimits_{i=1}^{N}\sum\nolimits_{j=1,j\neq
i}^{N}\big{|}\|\mathbf{g}_{i}-\mathbf{g}_{j}\|_{2}^{2}-\mathbf{D}_{i,j}^{2}\big{|},$
(20) $\displaystyle\text{s.t.}\
\eta_{1}\leq||\mathbf{g}_{i}-\mathbf{g}_{i+1}||_{2}\leq\eta_{2},\
\epsilon_{1}\leq||\mathbf{g}_{1}-\mathbf{g}_{N}||_{2}\leq\epsilon_{2},$

where $\eta_{1},\eta_{2}$ are defined in Eq. (1), $\epsilon_{1},\epsilon_{2}$
are defined in Eq. (2), $||\cdot||_{2}$ represents ${l}_{2}$ norm of a vector.

Then, in order to simplify the optimization process, we used the constraint
surrogate function $H(\cdot)$ (Eq. (13)) to turn the hard constraints in Eq.
(20) into soft constraints. Then we augment the learning objective and
reformulate an unconstrained optimization problem as follows,

$\displaystyle\underset{\mathbf{g}_{1},\cdots,\mathbf{g}_{N}\in\mathbb{R}^{3}}{\arg\min}\sum\nolimits_{i=1}^{N}\sum\nolimits_{j=1,j\neq
i}^{N}\big{|}\|\mathbf{g}_{i}-\mathbf{g}_{j}\|_{2}^{2}-\mathbf{D}_{i,j}^{2}\big{|}$
(21)
$\displaystyle+\lambda_{1}\sum_{i=1}^{N-1}H(||\mathbf{g}_{i}-\mathbf{g}_{i+1}||_{2};\eta_{1},\eta_{2},\eta_{3})+\lambda_{2}H(||\mathbf{g}_{1}-\mathbf{g}_{N}||_{2};\epsilon_{1},\epsilon_{2},\epsilon_{3}),$

where $H(\cdot)$ is the loss function defined in Eq. (13); hyperparameters
$\lambda_{1},\lambda_{2}\gg 1$ emphasize the importance of the constraints to
make sure that they are satisfied, hyperparameters $\eta_{3}$ and
$\epsilon_{3}$ are set to $\eta_{3}=\eta_{2}-\eta_{1}$,
$\epsilon_{3}=\epsilon_{2}-\epsilon_{1}$.

To optimize Eq. (21), we leverage random coordinate descent [24] to update
$\mathbf{g}_{1},\cdots,\mathbf{g}_{N}$ (i.e., randomly selecting a
$\mathbf{g}_{i}$ from $\\{\mathbf{g}_{1},\cdots,\mathbf{g}_{N}\\}$ to update
and fixing the remaining ones). $H(y;a,b,\delta)$ is convex w.r.t. $y$. Thus,
when we update each $\mathbf{g}_{i}$ ($i=1,\cdots,N$) in coordinate descent,
when fixing other $\mathbf{g}$ variables, Eq. (21) is convex w.r.t.
$\mathbf{g}_{i}$, which makes the optimization efficient. Worth mentioning
that the constrained 3D coordinates generation is not conducted during the
learning procedure and only during inference. Algorithm 1 summarizes the
essential steps.

##  3 Experiments

###  3.1 Experimental Setup

Dataset. (1) SAbDab (structural antibody database)
(http://opig.stats.ox.ac.uk/webapps/newsabdab/sabdab/) is a database
containing 5,494 antibodies structures [25]. (2) CoVAbDab (coronavirus
antibody database) [26] (http://opig.stats.ox.ac.uk/webapps/covabdab/)
documents 2,942 published binding antibodies to coronaviruses, the 3D
structures are not available. Data split strategies are given in Section E.

Metrics. (1) Protein Perplexity (PPL) measures protein likeness on amino acid
sequence [27]; (2) Root-Mean-Square Deviation (RMSD) measures the alignment
between generated 3D conformations and reference 3D conformation; (3) Validity
Rate (VR) is the percentage of valid (satisfying Eq.1, 2) CDR loops in all the
generated CDRs; (4) Diversity (Div) measures the diversity of generated CDR
loops on the amino acid sequence level. (5) SARS-CoV-2 Neutralization
Ability(used in [11]) is a well-trained machine learning model to evaluate the
neutralization ability of a CDR-H3 sequence to either SARS-CoV-1/2 viruses.
Detailed definitions are in Appendix.

Baselines. (1) Reference evaluates the results of CDR loops in the test set;
(2) LSTM (long short-term memory) [10, 9]; (3) GA (genetic algorithm) [8]; (4)
IR-GNN (Iterative refinement graph neural network) [11]; (5) AR-GNN (auto-
regressive graph neural network) [28, 29]. (6) AntibodyFlow \- constraint is
the variant of AntibodyFlow that does not consider validity constraint. LSTM,
IR-GNN and AR-GNN are deep generative models while GA is an evolutionary
algorithm. LSTM is unable to model 3D structures hence the results on RMSD and
VR are unavailable. All the approaches exhibit state-of-the-art performance on
antibody design [10, 9, 11]. More details are given in Appendix.

Table 1: CDR loop (H1, H2, H3) generation results on SAbDAb.

| Method |  PPL ($\downarrow$) |  RMSD ($\downarrow$) |  VR ($\uparrow$) |  Div ($\uparrow$)
---|---|---|---|---|---
H1 | Reference |  8.04$\pm$0.77 | 0.0 | 100.0% | 0.501
LSTM |  10.01$\pm$1.21 | - | - | 0.550
GA |  10.45$\pm$1.42 |  1.93$\pm$0.76 | 44.0% | 0.631
AR-GNN |  9.53$\pm$1.58 |  1.93$\pm$0.72 | 61.8% | 0.625
IR-GNN |  9.10$\pm$1.37 |  1.69$\pm$0.47 | 88.6% | 0.681
AntibodyFlow \- constraint |  9.32$\pm$1.21 |  1.92$\pm$0.75 | 71.4% | 0.635
AntibodyFlow | 8.68$\pm$0.92 | 1.28$\pm$0.45 | 96.3% | 0.645
H2 | Reference |  8.50$\pm$1.10 | 0.0 | 100.0% | 0.697
LSTM |  10.90$\pm$1.65 | - | - | 0.634
GA |  10.20$\pm$1.20 |  1.95$\pm$0.49 | 33.3% | 0.543
AR-GNN |  10.54$\pm$1.63 |  1.69$\pm$0.30 | 34.5% | 0.667
IR-GNN |  9.63$\pm$1.40 |  1.12$\pm$0.17 | 86.9% | 0.602
AntibodyFlow \- constraint |  9.96$\pm$1.89 |  1.70$\pm$0.30 | 49.5% | 0.672
AntibodyFlow | 9.23$\pm$1.32 | 0.99$\pm$0.11 | 92.8% | 0.674
H3 | Reference |  9.82$\pm$0.48 | 0.0 | 100.0% | 0.743
LSTM |  12.35$\pm$0.93 | - | - | 0.731
GA |  12.75$\pm$1.29 |  4.33$\pm$1.70 | 13.4% | 0.713
AR-GNN |  13.01$\pm$1.13 |  3.80$\pm$1.69 | 25.5% | 0.748
IR-GNN |  11.43$\pm$0.65 |  3.08$\pm$0.74 | 78.4% | 0.753
AntibodyFlow \- constraint |  10.96$\pm$0.80 |  3.21$\pm$0.98 | 59.1% | 0.769
AntibodyFlow | 10.81$\pm$0.47 | 2.50$\pm$0.50 | 91.0% | 0.780

Figure 2: Visualization of interpolation between two CDR H3 loops. For both
the 3D geometric structure and amino acid sequence, the changing trajectories
are smooth.  Table 2: SARS-CoV-2 neutralization optimization results (H3
sequence) on CoVAbDab.

Method |  PPL ($\downarrow$) |  Neutralization ($\uparrow$) |  Div ($\uparrow$)
---|---|---|---
Reference |  9.23$\pm$0.49 |  0.942$\pm$0.026 | 0.838
LSTM |  10.64$\pm$0.93 |  0.710$\pm$0.102 | 0.747
GA |  11.34$\pm$1.19 |  0.812$\pm$0.075 | 0.783
AR-GNN |  10.12$\pm$0.81 |  0.801$\pm$0.025 | 0.680
IR-GNN |  9.98$\pm$0.52 |  0.832$\pm$0.024 | 0.734
AntibodyFlow | 9.74$\pm$0.47 | 0.854$\pm$0.011 | 0.793

###  3.2 Results and Analysis

De novo design. We report the performance of de novo CDR H1, H2, and H3
generations in Table 1. For PPL and RMSD, we report the average value and
standard deviations of all the data points in the test set. We find that (i)
in all three tasks, our method achieves over 90% validity rate (96.3% on H1,
92.8% on H2, 91.0% on H3) and obtains the best performance among all the
methods in both sequences (PPL, perplexity) and geometric graph metrics (RMSD
& validity rate). Compared to IR-GNN (the best baseline), the relative
improvements are up to 5.4% in PPL (perplexity) (10.81 v.s. 11.43, H3), 24.3%
in RMSD (1.28 v.s. 1.69, H1), and 16.0% in VR (validity rate) (91.0% v.s.
78.4%, H3); (ii) among all the three tasks (H1, H2, H3), almost all the
methods get the highest perplexity, RMSD, and lower validity in the H3
generation task, thus validating the fact that CDR H3 loops have the highest
variability and most challenging to design [13, 11]; (iii) removing constraint
learning from AntibodyFlow (“AntibodyFlow \- constraint” v.s. AntibodyFlow)
significantly degrades the performance, especially in geometric graph metrics
(RMSD and validity rate (VR)). For example, in H3 loop generation (Table 1),
it leads to 1.4%, 22.1%, and 35.1% relative drops in PPL, RMSD, and VR,
respectively. It validates that learning constraint plays a crucial role in
geometric graph generation, especially generating a valid CDR loop; (iv)
AntibodyFlow achieves the highest diversity in H2 & H3 design tasks, and
second highest diversity in H1 generation, indicating our method is able to
generate diverse amino acid sequences.

Visualization of continuous latent variable. During the generation process,
the continuous latent variable $z$ controls the data space $x$ in the flow
model, i.e., $x=f^{-1}(z)$. Following [30, 19, 20], we visualize a linear
interpolation between two CDR loops to show their changing trajectory.
Specifically, given two CDR loops $x_{1}$ and $x_{2}$, we get the latent
variables $z_{1}=f(x_{1})$ and $z_{2}=f(x_{2})$. Then we get multiple $z$s
using linear interpolation between $z_{1}$ and $z_{2}$. We visualize these
$\tilde{z}$ via $\tilde{x}=f^{-1}(\tilde{z})$. We show the changing trajectory
of the generated CDR loops in Fig 2 and observe a smooth trajectory of both 3D
coordinates and amino acid sequence when varying the latent variable.

Optimization. The SARS-CoV-2 neutralization ability optimization results on
CoVAbDab data are reported in Table 2. 3D graphs are not available in CoVAbDab
data, hence we are unable to model the 3D graph. For PPL (perplexity) and
neutralization (SARS-CoV-2 neutralization ability), we report the average and
standard deviations. We find that AntibodyFlow obtains the lowest perplexity
(9.74 on average), highest neutralization ability (0.854 on average), and
highest diversity (0.793) among all the compared machine learning methods.

##  4 Conclusion

This paper proposed AntibodyFlow model to generate 3D antibody CDR loops.
Specifically, our method first builds the distance matrix of amino acids, and
then constructs the amino acid sequence conditioned on the distance matrix.
Also, to enhance the geometric validity of the generated CDR loops, we design
constraint learning and constrained 3D coordinate generation to impose the 3D
geometric constraints. Empirical studies validate AntibodyFlow’s superiority.

## Acknowledgment

This work is supported by the Fundamental Research Program of Shanxi Province
(No.20210302123025) and the Joint Funds of the Natural Science Foundation of
China (No.U21A20524).

## References

  * [1] Aaron L Nelson et al.  Development trends for therapeutic antibody fragments.  Nature biotechnology, 2009.
  * [2] Harini Narayanan et al.  Machine learning for biologics: opportunities for protein engineering, developability, and formulation.  Trends in pharmacological sciences, 2021.
  * [3] Robert M MacCallum et al.  Antibody-antigen interactions: contact analysis and binding site topography.  Journal of molecular biology, 1996.
  * [4] Inbal Sela-Culang, Vered Kunik, and Yanay Ofran.  The structural basis of antibody-antigen recognition.  Frontiers in immunology, 4:302, 2013.
  * [5] Titus Kretzschmar and Thomas Von Rüden.  Antibody discovery: phage display.  Current opinion in biotechnology, 2002.
  * [6] Eskil Söderlind et al.  Recombining germline-derived cdr sequences for creating diverse single-framework antibody libraries.  Nature biotechnology, 2000.
  * [7] RJ Pantazes and Costas D Maranas.  OptCDR: a general computational method for the design of antibody complementarity determining regions for targeted epitope binding.  Protein Engineering, Design & Selection, 2010.
  * [8] Xiaofeng Liu et al.  Computational design of an epitope-specific keap1 binding antibody using hotspot residues grafting and cdr loop swapping.  Scientific reports, 2017.
  * [9] Rahmad Akbar et al.  In silico proof of principle of machine learning-based antibody design at unconstrained scale.  BioRXiV, 2021.
  * [10] Koichiro Saka et al.  Antibody design using LSTM based deep generative model from phage display library for affinity maturation.  Scientific reports, 2021.
  * [11] Wengong Jin et al.  Iterative refinement graph neural network for antibody sequence-structure co-design.  ICLR, 2022.
  * [12] Victor Garcia Satorras, Emiel Hoogeboom, and Max Welling.  E(n) equivariant graph neural networks.  International Conference on Machine Learning, 2021.
  * [13] Cristian Regep et al.  The H3 loop of antibodies shows unique structural characteristics.  Proteins: Structure, Function, and Bioinformatics, 2017.
  * [14] Victor Garcia Satorras et al.  E(n) equivariant graph neural networks.  ICML, 2021.
  * [15] Yang Zhang.  I-tasser server for protein 3d structure prediction.  BMC bioinformatics, 2008.
  * [16] Robyn L Stanfield and Ian A Wilson.  Antibody structure.  Microbiology spectrum, 2(2):2–2, 2014.
  * [17] Jaroslaw Nowak et al.  Length-independent structural similarities enrich the antibody cdr canonical class model.  In MAbs. Taylor & Francis, 2016.
  * [18] Laurent Dinh et al.  NICE: Non-linear independent components estimation.  arXiv, 2014.
  * [19] Kaushalya Madhawa et al.  GraphNVP: An invertible flow model for generating molecular graphs.  arXiv, 2019.
  * [20] Chengxi Zang et al.  MoFlow: an invertible flow model for generating molecular graphs.  In KDD, 2020.
  * [21] Thomas N Kipf and Max Welling.  Semi-supervised classification with graph convolutional networks.  ICLR, 2016.
  * [22] Petar Veličković et al.  Graph attention networks.  ICLR, 2018.
  * [23] Congrui Yi and Jian Huang.  Semismooth newton coordinate descent algorithm for elastic-net penalized huber loss regression and quantile regression.  J. of Comp. and Graph. Stat., 2017.
  * [24] Ion Necoara et al.  A random coordinate descent algorithm for optimization problems with composite objective function and linear coupled constraints.  Computational Optimization and Applications, 2014.
  * [25] James Dunbar et al.  Sabdab: structural antibody database.  Nucleic acids research, 2014.
  * [26] Matthew IJ Raybould et al.  Cov-abdab: the coronavirus antibody database.  Bioinformatics, 2021.
  * [27] John Ingraham et al.  Generative models for graph-based protein design.  NeurIPS, 2019.
  * [28] Jiaxuan You et al.  Graphrnn: Generating realistic graphs with deep auto-regressive models.  In ICML, 2018.
  * [29] Shitong Luo et al.  A 3d generative model for structure-based drug design.  NeurIPS, 2021.
  * [30] Wengong Jin et al.  Junction tree variational autoencoder for molecular graph generation.  In ICML, 2018.
  * [31] Sam Sinai et al.  Variational auto-encoding of protein sequences.  NeurIPS MLCB workshop, 2017.
  * [32] Zak Costello and Hector Garcia Martin.  How to hallucinate functional proteins.  arXiv, 2019.
  * [33] Donatas Repecka et al.  Expanding functional protein sequence spaces using generative adversarial networks.  Nature Machine Intelligence, 2021.
  * [34] Mostafa Karimi et al.  De novo protein design for novel folds using guided conditional wasserstein generative adversarial networks.  Journal of Chemical Information and Modeling, 2020.
  * [35] Ge Liu et al.  Antibody complementarity determining region design using high-capacity machine learning.  Bioinformatics, 2020.
  * [36] Danilo Rezende et al.  Variational inference with normalizing flows.  In ICML, 2015.
  * [37] Diederik P Kingma and Max Welling.  Auto-encoding variational bayes.  arXiv, 2013.
  * [38] Christian Horvat and Jean-Pascal Pfister.  Denoising normalizing flow.  NeurIPS, 2021.
  * [39] Albert Pumarola et al.  C-flow: Conditional generative flow models for images and 3d point clouds.  In CVPR, 2020.
  * [40] Zineng Tang et al.  Continuous language generative flow.  In ACL, 2021.
  * [41] Chence Shi et al.  GraphAF: a flow-based autoregressive model for molecular graph generation.  In ICLR, 2020.
  * [42] Jose Juan Almagro Armenteros et al.  Language modelling for biological sequences–curated datasets and baselines.  BioRxiv, 2020.
  * [43] Wolfgang Kabsch.  A solution for the best rotation to relate two sets of vectors.  International Union of Crystallography, 1976.
  * [44] Zhiheng Huang, Wei Xu, and Kai Yu.  Bidirectional LSTM-CRF models for sequence tagging.  arXiv, 2015.

###### Contents

  1. 1 Introduction
  2. 2 Method
    1. 2.1 Problem Formulation
    2. 2.2 Background and Overview
    3. 2.3 Distance Matrix Flow
    4. 2.4 Conditional Amino Acid Flow
    5. 2.5 Constraint Learning
    6. 2.6 Constrained 3D Coordinates Generation
  3. 3 Experiments
    1. 3.1 Experimental Setup
    2. 3.2 Results and Analysis
  4. 4 Conclusion
  5. A Overview of Appendix
  6. B Related Work
  7. C Proof of Corollary 1
  8. D Learning and Inference
  9. E Data and Preprocessing
  10. F Evaluation Metrics
  11. G Baseline Methods
  12. H Implementation Details
  13. I Sensitivity Analysis
  14. J Case Study

## Appendix

##  Appendix A Overview of Appendix

The supplementary materials are organized as follows. Table 3 lists all the
important mathematical notations. Section B reviews the existing literature.
Section C proves Corollary 1. We describe the learning and inference process
in Section D. Then we elaborate on the experimental setup, including data and
preprocessing in Section E, evaluation metrics in Section F, and baseline
methods in Section G. Implementation details are given in Section H. Section I
conducts sensitivity analysis on some key hyperparameters. Section J shows a
case study. We also upload the code repository in the supplementary
materials222Code repository is available in the supplementary materials. .

Figure 3: The whole framework of AntibodyFlow. (Forward) Flow model $f$
incorporates (i) ${\mathbf{z}}_{\mathbf{D}}=f_{\mathbf{D}}(\mathbf{D})$ (Sec.
2.3) and (ii)
${\mathbf{z}}_{\mathbf{S}|\mathbf{D}}=f_{\mathbf{S}|\mathbf{D}}(\mathbf{S}|\mathbf{D})$
(Sec. 2.4); Inverse flow model $f^{-1}$ incorporates (i)
$D=f_{\mathbf{D}}^{-1}({\mathbf{z}}_{\mathbf{D}})$ (Sec 2.3) and (ii)
$S=f_{\mathbf{S}|\mathbf{D}}^{-1}({\mathbf{z}}_{\mathbf{S}|\mathbf{D}};\mathbf{z}_{\mathbf{D}})$
(Sec 2.4); “Constraint learning” (Sec. 2.5) minimizes constraint loss (Eq. 17)
to encourage flow model to learn these constraints. “Constrained 3D
coordinates generation” (Sec 2.6) generates 3D coordinates based on distance
matrix and validity constraints of CDR loops defined in Eq. (1) and (2).

##  Appendix B Related Work

(I) de novo antibody/protein design is about constructing antibodies/proteins
from scratch with desirable properties. For protein design, variational
autoencoder (VAE) based methods [31, 32] and generative adversarial network
(GAN) based methods [33, 34] were proposed to generate novel amino acid
sequences. For antibody design, [7] combines CDR canonical structures and
iteratively redesigns their positions; [35] designed a neural network ensemble
to backpropagate the gradient to update the amino acid sequence in continuous
space; [10, 9] leverage Long Short Term Memory (LSTM) to generate amino acid
sequences; [11] proposed an iterative refinement graph neural network to
jointly design 3D graph structure and amino acid sequence. However, almost all
of these methods incrementally grow the sequence by adding a single amino acid
at a time. In contrast, our method generates the whole CDR loop together so
that the global information (geometric constraints) is considered.

(II) Flow model (a.k.a. normalizing flow model) is a special category of deep
generative models that designs invertible mapping so that the input data can
be exactly reconstructed [36]. The advantage of a flow model is its ability to
reconstruct input data exactly and directly optimize the likelihood, instead
of optimizing a surrogate, e.g., VAE optimizes the lower bound of likelihood
[37]. It is successfully applied to image generation/impainting [38, 39],
machine translation [40], and small molecule generation [41, 19, 20, 14]. The
conditional flow model is a well-known variant of the flow model and is also
widely applied [39, 19, 20]. In this paper, we address some challenges of the
flow model in 3D graph generation. Besides, [39] proposes conditional flow on
computer vision.

(III) 3D graph generative model have been proposed for drug discovery
recently, e.g., [11] (mentioned in I) and [29]. [29] leveraged autoregressive
3D graph neural network to design small-molecule conformation conditioned on
the pocket conformation of the target protein. However, these methods focus on
generating a single component (e.g., an amino acid or an atom) at a time based
on conditional distribution without global constraints. In contrast, our
method generates the 3D graph at once and imposes global constraints that
encode the geometric shape of CDR loops.

Table 3: Mathematical notations.

Notations | Explanations
---|---
$\mathbf{x}=(\mathbf{D},\mathbf{S})$ | 3D CDR loop
$N$ | number of amino acids in CDR loop
$V,|V|=20$ | vocabulary of amino acids
$\mathbf{G}=[\mathbf{g}_{1},\cdots,\mathbf{g}_{N}]^{\top}$ |  3D coordinates of amino acids, $\mathbf{G}\in\mathbb{R}^{N\times 3}$
$\mathbf{S}\in\\{0,1\\}^{N\times|V|}$ | amino acid matrix
$\mathbf{D}\in\mathbb{R}_{+}^{N\times N}$ | pairwise distance matrix
$f$ (i.e., $f_{\theta}$), |  flow model, learnable parameters $\theta$
$\frac{\partial f_{\theta}}{\partial\mathbf{x}}$, i.e., $\frac{\partial{\mathbf{z}}}{\partial\mathbf{x}}$ |  Jacobian of $f$
$\det$ | determinant
$\eta_{1},\eta_{2},\epsilon_{1},\epsilon_{2}$ |  validity constraint, Eq. (1), (2)
$H(y;a,b,\delta)$ |  constraint surrogate function, Eq. (13)
$h(\mathbf{D})$ |  constraint loss function of $\mathbf{D}$, Eq. (17)
$\xi_{1},\xi_{2},\xi_{3}$ / $\lambda_{1},\lambda_{2}$ |  hyperparameter in Eq. (17) / (21)

##  Appendix C Proof of Corollary 1

Corollary 1. Suppose $\mathbf{G}^{\prime}\in\mathbb{R}^{N\times 3}$ are the 3D
coordinates of the graph by rotating and translating the original graph whose
coordinates are $\mathbf{G}\in\mathbb{R}^{N\times 3}$,
$\mathbf{D}\in\mathbb{R}^{N\times N}$ and
$\mathbf{D}^{\prime}\in\mathbb{R}^{N\times N}$ are the distance matrix of
$\mathbf{G}$ and $\mathbf{G}^{\prime}$, respectively, then we have
$\mathbf{D}=\mathbf{D}^{\prime}$.

###### Proof.

Suppose $\mathbf{D}$ and $\mathbf{D}^{\prime}$ are the distance matrix before
and after transformation (rotation/translation), respectively.

First, regarding rotation, without loss of generalization, we suppose the
rotation matrix is $\mathbf{Q}\in\mathbb{R}^{3\times 3}$, which is an
orthogonal matrix satisfying
$\mathbf{Q}^{\top}\mathbf{Q}=\mathbf{Q}\mathbf{Q}^{\top}=\mathbf{I}_{3}$,
$\mathbf{I}_{3}\in\mathbb{R}^{3\times 3}$ is identity matrix. Rotating 3D
graph whose coordinates are
$\mathbf{G}=[\mathbf{g}_{1},\cdots,\mathbf{g}_{N}]$ results in a new 3D graph
whose coordinates are
$\mathbf{G}^{\prime}=[\mathbf{Q}\mathbf{g}_{1},\cdots,\mathbf{Q}\mathbf{g}_{N}]$.
Then based on the definition of the distance matrix, we have

$\displaystyle\mathbf{D}^{\prime}_{i,j}$
$\displaystyle=||\mathbf{Q}\mathbf{g}_{i}-\mathbf{Q}\mathbf{g}_{i}||_{2}$ (22)
$\displaystyle=||\mathbf{Q}(\mathbf{g}_{i}-\mathbf{g}_{j})||_{2}$
$\displaystyle=\sqrt{(\mathbf{g}_{i}-\mathbf{g}_{j})^{\top}\mathbf{Q}^{\top}\mathbf{Q}(\mathbf{g}_{i}-\mathbf{g}_{j})}$
$\displaystyle=\sqrt{(\mathbf{g}_{i}-\mathbf{g}_{j})^{\top}(\mathbf{g}_{i}-\mathbf{g}_{j})}$
$\displaystyle=||\mathbf{g}_{i}-\mathbf{g}_{j}||_{2}$
$\displaystyle=\mathbf{D}_{i,j}.$

Secondly, translating 3D graph (whose coordinates are
$\mathbf{G}=[\mathbf{g}_{1},\cdots,\mathbf{g}_{N}]$) by
$\mathbf{v}\in\mathbb{R}^{3}$ results in a new 3D graph whose coordinates are
$\mathbf{G}^{\prime}=[\mathbf{g}_{1}+\mathbf{v},\cdots,\mathbf{g}_{N}+\mathbf{v}]$.
We have

$\mathbf{D}^{\prime}_{i,j}=||(\mathbf{g}_{i}+\mathbf{v})-(\mathbf{g}_{j}+\mathbf{v})||_{2}=||\mathbf{g}_{i}-\mathbf{g}_{j}||_{2}=\mathbf{D}_{i,j}.$
(23)

Third, the joint transformation of rotation and translation can be represented
as

$\mathbf{G}^{\prime}=[\mathbf{Q}\mathbf{g}_{1}+\mathbf{v},\cdots,\mathbf{Q}\mathbf{g}_{N}+\mathbf{v}].$

Then based on definition of distance matrix and combining Equation (22) and
(23), we have

$\displaystyle\mathbf{D}^{\prime}_{i,j}$
$\displaystyle=||(\mathbf{Q}\mathbf{g}_{i}+\mathbf{v})-(\mathbf{Q}\mathbf{g}_{i}+\mathbf{v})||_{2}$
(24) $\displaystyle=||\mathbf{Q}\mathbf{g}_{i}-\mathbf{Q}\mathbf{g}_{i}||_{2}$
$\displaystyle=\mathbf{D}_{i,j}.$

In sum, the distance matrix $\mathbf{D}$ is invariant w.r.t. translation and
rotation on the 3D graph.

∎

##  Appendix D Learning and Inference

In this section, we describe the learning and inference process. The key steps
of our method are summarized in Algorithm 1.

During the learning process, we alternatively conduct maximum likelihood
learning and constraint learning. Given the data sample $\mathbf{x}$, maximum
likelihood learning (MLE) maximizes the exact log-likelihood (Equation. 3) of
input data w.r.t. parameters $\theta$,

$\displaystyle\mathcal{L}_{\text{MLE}}(\theta)=\log
p_{{\mathbf{x}}}({\mathbf{x}})$ (25) $\displaystyle=$ $\displaystyle\log
p_{\mathbf{z}}({\mathbf{z}})+\log\big{(}\big{|}\det\frac{\partial{\mathbf{z}}}{\partial{\mathbf{x}}}\big{|}\big{)}$
$\displaystyle=$ $\displaystyle\log
p_{{\mathbf{z}}_{\mathbf{D}}}({\mathbf{z}}_{\mathbf{D}})+\log
p_{{\mathbf{z}}_{\mathbf{S}|\mathbf{D}}}({\mathbf{z}}_{\mathbf{S}|\mathbf{D}})$
$\displaystyle+\log\bigg{(}\bigg{|}\det\frac{\partial
f_{\mathbf{D}}}{\partial\mathbf{D}}\bigg{|}\bigg{)}+\log\bigg{(}\bigg{|}\det\frac{\partial
f_{\mathbf{S}|\mathbf{D}}}{\partial\mathbf{S}}\bigg{|}\bigg{)},$

where latent variable $\mathbf{z}$ satisfies
$\mathbf{z}=f_{\theta}(\mathbf{x})$, the first two term
$p_{{\mathbf{z}}_{\mathbf{D}}}(\cdot)$,
$p_{{\mathbf{z}}_{\mathbf{S}|\mathbf{D}}}(\cdot)$ follow isotropic Gaussian
distribution, the third term $\log\big{(}\big{|}\det\frac{\partial
f_{\mathbf{D}}}{\partial\mathbf{D}}\big{|}\big{)}$ is given in Section 2.3;
the fourth term $\log\big{(}\big{|}\det\frac{\partial
f_{\mathbf{S}|\mathbf{D}}}{\partial\mathbf{S}}\big{|}\big{)}$ is defined in
Section 2.4.

During inference process, we sample latent variable
${\mathbf{z}}=({\mathbf{z}}_{\mathbf{D}},{\mathbf{z}}_{\mathbf{S}|\mathbf{D}})$
from isotropic Gaussian distribution, generate the distance matrix
$\mathbf{D}$, and transform $\mathbf{D}$ into 3D coordinates $\mathbf{G}$
following Section 2.3. Then we generate $\mathbf{S}$ conditioned on
$\mathbf{D}$. For each row in $\mathbf{S}\in\mathbb{R}^{N\times|V|}$, we use
argmax operation to the $|V|$-dimensional row vector to obtain the
corresponding amino acid’s category. The whole pipeline is illustrated in
Figure 3.

##  Appendix E Data and Preprocessing

This section describes the datasets and the preprocessing procedure.

(1) SabDab (the structural antibody database)333It is publicly available
http://opig.stats.ox.ac.uk/webapps/newsabdab/sabdab/.  is a database
containing all the antibody structures available in the PDB (Protein Data
Bank), annotated and presented in a consistent fashion [25]. It has 5,494
antibodies. One antibody may contains multiple CDR loops (both amino acid
sequence and 3D structure). After filtering the antibodies whose CDR
structures are not available, we get 8,381 data points. For H1 loops, the
average length of CDR loops (i.e., number of amino acids) is 7.2, the standard
deviation is 0.6; for H2 loops, the average length is 5.9, with a 0.7 standard
deviation; for H3 loops, the average is 12.3, and the standard deviation is
4.4. The data points are randomly split into training, validation, and test
sets with a ratio of 8:1:1. The data is saved in the format of a PDB (Protein
Data Bank) file, where the coordinates of all the atoms are available. We use
the coordinate of $\alpha$-carbon atom to denote the coordinate of amino acid,
following [27, 11]. The unit of distance in PDB file is $\mathring{A}$,
$1\mathring{A}=10^{-10}m=0.1nm$.

(2) CoV-AbDab (the coronavirus antibody database) [26]444It is publicly
available at http://opig.stats.ox.ac.uk/webapps/covabdab/ documents all
published/patented binding antibodies and nanobodies to coronaviruses,
including SARS-CoV2, SARS-CoV1, and MERS-CoV. It contains 2,942 antibodies and
494 nanobodies. The 3D structures are not available. In this paper, we
restrict our attention on CDR H3 sequence of antibodies, and ignore
nanobodies. We use half of the data (1.5K) to train machine learning model to
evaluate SARS-CoV-2 neutralization ability (as mentioned in Section F), it is
a binary classification task, the train/validation/test split ratio is 8:1:1.
The remaining data are used to train H3 loops generation models, where
train/validation split ratio is also 9:1. Test data is not required because
the generated H3 loops are evaluated using the trained model that evaluates
SARS-CoV-2 neutralization ability.

##  Appendix F Evaluation Metrics

In this section, we describe the metrics for evaluating the performance of all
the methods.

  * •

(1) Protein Perplexity (PPL). Protein perplexity is the exponential of the
average log-likelihood and measures how well a probability model predicts
amino acid sequence [27, 42]. We train an LSTM model as language model. Lower
perplexity indicates better performance.

  * •

(2) Root-Mean-Square Deviation (RMSD) measures the alignment between tested
conformations $R\in\mathbb{R}^{N\times 3}$ and reference conformation
${R}^{r}\in\mathbb{R}^{N\times 3}$, $N$ is the number of points in the
conformation, defined as

$\text{RMSD}(R,\hat{R})=\big{(}\frac{1}{N}\sum_{i=1}^{N}||R_{i}-\hat{R}_{i}||_{2}^{2}\big{)}^{\frac{1}{2}},$

where $R_{i}$ is the $i$-th row of conformation $R$. The conformation
$\hat{R}$ is obtained by an alignment function $\hat{R}=A(R,R^{r})$, which
rotates and translates the reference conformation $R^{r}$ to have the smallest
distance to the generated R according to the RMSD metrics, which is calculated
by the Kabsch algorithm [43]. Lower RMSD indicates better performance.

  * •

(3) Validity Rate (VR) is the percentage of valid CDR loops in all the
generated CDRs. A CDR loop is valid if it satisfies two criteria in Equation
(1) and (2). Higher validity is more desirable.

  * •

(4) Similarity is evaluated on amino acid sequence level. The similarity
between two sequences $s_{1}$ and $s_{2}$, and is defined as

$\text{sim}(s_{1},s_{2})=\frac{|\text{LCS}(s_{1},s_{2})|}{\max\\{|s_{1}|,|s_{2}|\\}},$
(26)

where $|s|$ denotes the length of sequence $s$, $\text{LCS}(s_{1},s_{2})$
denotes the longest common substring between $s_{1}$ and $s_{2}$, computed by
dynamic programming. The similarity value ranges from 0 to 1, a higher value
indicates more similarity. The distance is defined as one minus similarity.
The similarity is not used directly in our experiment. It is used for
evaluating diversity.

  * •

(5) Diversity (Div) is defined as the average pairwise distance between the
CDR loops. Specifically, suppose $\mathcal{S}$ is the set of all the generated
CDR loops, we have

$\text{diversity}(\mathcal{S})=1-\frac{1}{|\mathcal{S}|(|\mathcal{S}|-1)}\sum_{s_{1},s_{2}\in\mathcal{S},s_{1}\neq
s_{2}}\text{sim}(s_{1},s_{2}),$

where $\text{sim}(\cdot,\cdot)$ is the similarity function defined above.
Diversity ranges from 0 to 1, higher diversity is more desirable.

  * •

(6) SARS-CoV Neutralization Ability. Following [11], we train a machine
learning model to predict the neutralization ability of a CDR-H3 sequence to
either SARS-CoV-1 or SARS-CoV-2 viruses. The predictor takes as CDR-H3
sequence of an antibody and outputs a neutralization probability for the SARS-
CoV-1 and SARS-CoV-2 viruses. Each residue is embedded into a 100-dimensional
vector. Then the embeddings are fed to a three-layer bi-LSTM (bidirectional
Long Short Term Memory) [44], followed by average-pooling and a two-layer MLP
(multi-layer perceptron). The hidden dimension of bi-LSTM is 100. The activate
functions in the hidden layers and output layer are ReLU and sigmoid,
respectively. The prediction is very accurate, achieving over 0.99 ROC-AUC and
over 0.95 F1 score. The neutralization predictor yields neutralization
probability (ranging from 0 to 1), higher probability is desirable.

##  Appendix G Baseline Methods

In this section, we describe the baseline methods and briefly introduce the
implementation details. All these baselines are state-of-the-art approaches in
antibody design. We use the default setup following the original paper.

  1. 1.

Reference evaluates the metrics of real CDR loops in the test set.

  2. 2.

LSTM [10, 9] leverages LSTM (long short-term memory) to generate amino acid
sequence and does not model 3D structure. Following the original paper, it is
trained by the Adam optimizer and has a learning rate of 0.01. Dropout rates
were chosen from 0.1 and 0.2 to regularize all layers The hidden dimension of
LSTM is set to 256.

  3. 3.

GA (genetic algorithm) starts from a population of CDR loops. Given the data
in the current population, to generate offspring for the next generation, we
leverage the following two operations: (i) crossover (switching subsequence
between two loops) and (ii) mutation (randomly changing one amino acid) to
produce new candidate in each iteration [8]. The population size of GA is set
to 100, and the number of generations is set to 1K, which is sufficient for GA
to converge. In each generation, when producing the offspring, we generate 1K
offsprings, where half of them are produced through crossover operation, and
half are produced through mutation operation. We use perplexity as the oracle
to select the most promising offspring for the next generation.

  4. 4.

IR-GNN (Iterative refinement graph neural network) [11] first unravels amino
acid sequence and iteratively refine its predicted global structure. The
inferred structure in turn guides subsequent amino acid choices. Following
their original paper, for the GNN, both its structure and sequence MPN have
four message passing layers, with a hidden dimension of 256 and block size b =
4. It is trained by the Adam optimizer with a dropout of 0.1 and a learning
rate of 0.001.

  5. 5.

AR-GNN (auto-regressive graph neural network) [28, 29] was proposed for
molecule generation and was adapted for antibody CDR design [11]. The neural
architecture of AR-GNN follows the setups of IR-GNN. Different from IR-GNN
that uses iterative refinement mechanism, AR-GNN generates the amino acid
sequence and 3D structural graph autoregressively.

  6. 6.

AntibodyFlow \- constraint is the variant of AntibodyFlow that does not
consider validity constraint (Section 2.5 and 2.6). The other setups are the
same as AntibodyFlow (Section H).

LSTM, IR-GNN and AR-GNN fall into the category of deep generative models while
GA are belonging to traditional heuristic searching method. [14] is an ODE
type of flow model. However, the training process is computationally expensive
since the same forward operation has to be done multiple times sequentially in
order to solve the ODE equation. It also exhibited some instabilities [14] and
fails in our experiment, so it is not included as baseline. All the approaches
exhibit state-of-the-art performance on antibody design [10, 9, 11].

##  Appendix H Implementation Details

In this section, we elaborate on the implementation details for future
reproduction, including network architecture, hyperparameter setting, etc. It
is implemented using Pytorch 1.7.1, and Python 3.7.1. The model is trained on
NVIDIA Pascal Titan X GPU. We also upload the code repository in the
supplementary materials.

We use 5 affine coupling layers for the distance matrix flow model and 10
affine coupling layers for the conditional amino acid flow model. We adopt two
3 $\times$ 3 convolution layers with 128 dimensions in the hidden state of
every coupling layer. For each graph coupling layer, we set one wGNN layer
(Sec 2.4) with 64 dimensions followed by a two-layer multilayer perceptron
with 128 & 64 hidden dimensions. $\xi$ in Equation (12) is set to 0.3. The
model is trained in an end-to-end manner, where we use Adam optimizer with
initial learning rate $1e^{-3}$, the batch size is set to 64, and the number
of training epochs is set to 200. The model is trained on a NVIDIA Pascal
Titan X GPU, the training process takes around 3 hours for each of the three
tasks on SAbDab.

In constraint learning (Section 2.5), the hyperparameters
$\xi_{1},\xi_{2},\xi_{3}$ in Equation (17) are set to
$\xi_{1}=10,\xi_{2}=50,\xi_{3}=1$. We alternatively optimize log-likelihood in
Equation. (3) and constraint loss in Equation. (17). Constraint learning is
only conducted during the learning procedure and is not utilized in inference.

Figure 4: An illustration of H function (Eq. 13), where $a$=$-1$, $b$=$1$,
$\delta$=$1$. It has equal values and slopes of the different sections at the
four connection points ($-2,-1,1,2$) so it is differentiable everywhere.

Figure 5: Sensitivity Analysis on $\lambda_{1}$ and $\lambda_{2}$. We find the
combination $\lambda_{1}=50$ and $\lambda_{2}=100$ works best empirically.

In 3D coordinates generation (Equation. 21 in Section 2.6), to make sure that
the constraints are satisfied, we set $\lambda_{1}=50$ and $\lambda_{2}=100$.
We conduct sensitivity analysis about $\lambda_{1}$ and $\lambda_{2}$ in
Section I to justify their choices. Hyperparameters $\eta_{3}$ and
$\epsilon_{3}$ are set to $\eta_{3}=\eta_{2}-\eta_{1}$,
$\epsilon_{3}=\epsilon_{2}-\epsilon_{1}$. We use SGD (stochastic gradient
descent) as optimizer when optimizing Equation (21). Worth to mention that the
constrained 3D coordinates generation is not conducted during the learning
procedure and is only conducted during inference.

The setup of $\eta_{1},\eta_{2},\epsilon_{1},\epsilon_{2}$ is based on domain
knowledge [16, 17] and empirical validation. Concretely, for H1 loop,
$\eta_{1}=3.76$, $\eta_{2}=3.84$, $\epsilon_{1}=11.4,\epsilon_{2}=13.1$; for
H2 loop, $\eta_{1}=3.76$, $\eta_{2}=3.87$,
$\epsilon_{1}=5.0,\epsilon_{2}=5.9$; for H3 loop, $\eta_{1}=3.71$,
$\eta_{2}=3.88$, $\epsilon_{1}=6.5,\epsilon_{2}=8.5$. The units of
$\epsilon_{1},\epsilon_{2},\eta_{1},\eta_{2}$ are Angstrom $\mathring{A}$
($10^{-10}$ m, i.e., 0.1 nm). We find that the bond lengths between amino
acids are relatively fixed across H1, H2, and H3 loops, whereas the lengths of
open loops vary greatly. The hyperparameters $\epsilon_{3},\eta_{3}$ are set
to $\epsilon_{3}=\epsilon_{2}-\epsilon_{1}$ and $\eta_{3}=\eta_{2}-\eta_{1}$.

On each generation task (H1, H2, H3), we generate 10K CDR loops and evaluate
their VR (validity rate), PPL (perplexity), and Div (Diversity). When
evaluating geometric graph metrics RMSD, for each tested CDR loop in the test
set, we select the CDR loop that has the smallest RMSD with the test CDR loops
from the set of all generated CDR loops and evaluate the RMSD.

Figure 6: Case study of CDR H3 structure of Reference, AntibodyFlow, IR-GNN
(the best baseline). The H3 loop is from the antibody whose PDB ID is “3s34”.
We found that (i) on amino acid sequence level, the similarity (defined in
Equation 26) between “VTDAFMI” (generated by AntibodyFlow) and “VTDAFDI” (the
reference) is 0.86 ($\frac{6}{7}$), whereas the similarity between “VTYAFRG”
(IR-GNN) and “VTDAFDI” (the reference) is 0.57. Also, AntibodyFlow got lower
RMSD and perplexity than IR-GNN.

##  Appendix I Sensitivity Analysis

This section conducts sensitive analysis on some critical hyperparameters. In
3D coordinates generation (Equation. 21 in Section 2.6), we use $\lambda_{1}$
and $\lambda_{2}$ to incorporate soft constraints to make sure that the hard
constraints (defined in Equation 1 and 2) are satisfied.

We conduct sensitive analysis about $\lambda_{1}$ and $\lambda_{2}$ on cdr H1
loop design tasks. Specifically, regarding $\lambda_{2}$ $\lambda_{2}$, we fix
one and vary the other using grid search, then record the results in terms of
RMSD (root mean square deviation) and validity rate. In Figure 5, we plot the
RMSD and validity as a function of $\lambda_{1}$ (fixing $\lambda_{2}$) and
$\lambda_{2}$ (fixing $\lambda_{1}$), respectively. We find the combination
$\lambda_{1}=50$ and $\lambda_{2}=100$ works best empirically. We apply the
similar strategy to select other hyperparameters, e.g.,
$\xi_{1},\xi_{2},\xi_{3}$ in Equation 17.

##  Appendix J Case Study

This section shows a case study in Figure 6 and analyzes the results.
Specifically, we show CDR H3 structure of Reference, AntibodyFlow, IR-GNN (the
best baseline). The H3 loop is from the antibody whose PDB ID is “3s34”. We
found that (i) on amino acid sequence level, the similarity (defined in
Equation 26) between “VTDAFMI” (generated by AntibodyFlow) and “VTDAFDI” (the
reference) is 0.86 ($\frac{6}{7}$), whereas the similarity between “VTYAFRG”
(IR-GNN) and “VTDAFDI” (the reference) is 0.57. Also, AntibodyFlow got lower
RMSD and perplexity than IR-GNN. The case study validates AntibodyFlow’s
superiority.