#  Advancing Model Pruning via Bi-level Optimization

David S. Hippocampus
Department of Computer Science
Cranberry-Lemon University
Pittsburgh, PA 15213
<EMAIL_ADDRESS>
Use footnote for providing further information about author (webpage,
alternative address)— _not_ for acknowledging funding agencies.

###### Abstract

The deployment constraints in practical applications necessitate the pruning
of large scale deep learning models, i.e., promoting their weight sparsity. As
illustrated by the Lottery Ticket Hypothesis (LTH), pruning also has the
potential of improving their generalization ability. At the core of LTH,
iterative magnitude pruning (IMP) is the predominant pruning method to
successfully find “winning tickets”, i.e., pruned sparse models with better
generalization than the original dense models. Yet, the computation cost of
IMP prohibitively grows as the targeted pruning ratio increases. To reduce the
computation overhead, various efficient ‘one-shot’ pruning methods have been
developed but these schemes are usually unable to find winning tickets as good
as IMP. This raises the question of how to close the gap between pruning
accuracy and pruning efficiency? To tackle it, we pursue the algorithmic
advancement of model pruning. Specifically, we formulate the pruning problem
from a fresh and novel viewpoint, bi-level optimization (BLO). We show that
the BLO interpretation provides a technically-grounded optimization base for
an efficient implementation of the pruning-retraining learning paradigm used
in IMP. We also show that the proposed bi-level optimization-enabled pruning
method (termed BiP) is a special class of BLO problems with a bi-linear
problem structure. By leveraging such bi-linearity, we theoretically show that
BiP can be solved as easily as first-order optimization, thus owning the
computation efficiency. Through extensive experiments on both structured and
unstructured pruning with 5 model architectures and 4 data sets, we
demonstrate that BiP can find better winning tickets than IMP in most cases,
and is computationally as efficient as the one-shot pruning schemes,
demonstrating 2-7$\times$ speedup over IMP for the same level of model
accuracy and sparsity.

##  1 Introduction (1 page)

Deep neural networks (DNNs) have dramatically revolutionized the state-of-the-
art performance in various fields [krizhevsky2012imagenet, brown2020language].
As a key to increase the generalization of the DNNs [bartlett2021deep], the
over-parameterized structures also bring its own problem. The millions or even
billions of parameters not only increase computational costs during inference,
but also raise great challenges for deployment on the resource-limited devices
[han2015deep]. Thus, the model pruning techniques has become the research hot
spot recent years to remove redundant parameters and reduce the model size
[han2015deep, blalock2020state].

Among numerous studies on model pruning algorithms [lecun1989optimal,
ren2018admmnn, ramanujan2020sHiddenSubnetwork], the heuristics-based Iterative
Magnitude Pruning (IMP) is the current dominant approach to achieve model
sparsity [frankle2018lottery, Renda2020Comparing, frankle2020linear,
blalock2020state, chen2020lottery, chen2021lottery], while the Lottery Ticket
Hypothesis (LTH) serves as one of the most well-known approaches to attain
sparse sub-networks without suffering performance loss. LTH hypothesizes the
existence of a sub-network (the so-called ‘winning ticket’), when trained from
the scratch, can match the performance of the full model. The core idea of LTH
lies in repetitively calling the pruning-retraining cycle and prunes a small
ratio of the remaining parameters at each round. However, the achieved winning
ticket is often at the cost of even more computations than training a dense
model from scratch, which hinders its scalability towards large dataset and
large model architectures. 11todo: 1Is LTH an algorithm or a hypothesis that
motivates algorithms like IMP? To improve the efficiency of LTH, numerous
heuristics-based one-shot pruning methods [ma2021sanity, frankle2018lottery,
wang2020pick, lee2018snip, tanaka2020pruning, alizadeh2022prospect], _e.g._ ,
one-shot magnitude pruning (OMP), have been proposed, which directly prunes
the model to the target sparsity and achieves significant speedup compared to
IMP. However, all current the one-shot methods encountered drastic performance
drop and the odds of finding the winning tickets are also small. 22todo:
2substantiate with citations or support otherwise

Another research direction is the optimization-based pruning methods
[sehwag2020hydra, ramanujan2020sHiddenSubnetwork, lee2020layer,
csordas2020neural, ren2018admmnn], which still follow the pruning-retraining
paradigm and adopt sparsity regularization [liu2019admm, wang2020neural] or
parameterized masking [csordas2020neural, ramanujan2020sHiddenSubnetwork,
sehwag2020hydra] to promote sparsity. However, current optimization-based
pruning still suffers from a non-negligible performance gap from IMP and can
hardly find the winning tickets. 33todo: 3substantiate with citations or
support otherwise This motivates us to rethink the fundamentals of this
research direction: can we advance the optimization basis of the optimization-
based pruning and achieve comparable or even better performance than IMP?

Figure 1:  A performance snapshot of the proposed BiP method vs. the IMP
baseline and the original dense model across three pruned ResNet models
(ResNet-20, ResNet-56, and ResNet-18) with 74% sparsity on CIFAR-10. The
marker size indicates the model size (i.e., the remaining model parameters).
And the uni-color region corresponds to the same model type used by different
pruning methods.

To tackle the above challenge, we revisit the pruning problem from a fresh and
novel viewpoint of bi-level optimization (BLO). In the new formulation, we
show that BLO provides a technically-grounded optimization basis for an
efficient implementation of the pruning-retraining learning paradigm.

To the best of our knowledge, this is the first work making a rigorous
connection between pruning and BLO. Technically, we propose a novel bi-level
optimization-enabled pruning method (termed BiP). We further show BiP takes
advantage of the bi-linearity of the pruning variables and avoids the
computationally intensive nature of the common BLO methods. More specifically,
BiP can be solved as easily as the first-order alternating optimization
problem. Practically, extensive experiments demonstrate the superiority of the
proposed BiP in terms of accuracy, sparsity, and efficiency. BiP achieves the
best winning tickets in all settings with only OMP-level time consumption. As
Fig. 1 shows, BiP achieves winning tickets with 5$\times$ speedup compared to
IMP on CIFAR-10 with ResNet-18 in the sparsity of 74%. 44todo: 4MINOR: in fig
1 make axis labels, legend and arch name larger

Our main contributions are listed as the following aspects:

$\bullet$ (Formulation-wise) We reformulate the model pruning problem as a bi-
level optimization (BLO) problem. We show that the new formulation not only
advances the optimization foundation of the pruning problem, but also provides
a natural interpretation of the widely-accepted pruning-retraining paradigm.
More importantly, the new formulation disentangles model pruning and
retraining variables and provides more flexibility to design mismatched
objectives at the upper and lower optimization levels.

$\bullet$ (Algorithmic-wise), we propose a new pruning algorithm BiP based on
our new BLO formulation and the implicit gradient theory. Different from the
ordinary computationally intensive BLO methods, we theoretically show that BiP
can be solved as easily as the first-order optimization by taking advantage of
the bi-linear nature of the pruning variables.

$\bullet$ (Experiment-wise) We conduct extensive experiments across 4 datasets
(CIFAR-10/100, Tiny-ImageNet and ImageNet), 5 model architectures, and 2
pruning settings (unstructured and structured pruning). We show that BiP
combines the high accuracy of the predominant iterative magnitude pruning
(IMP) and the high efficiency of the one-shot pruning schemes. Specifically,
(i) BiP achieves higher test accuracy than IMP and finds the best winning
tickets in all settings. (ii) BiP enjoys high pruning efficiency (comparable
to one-shot pruning schemes), and achieves $2\sim 7\times$ speedup over IMP
for the same level of model accuracy and sparsity.

##  2 Related Work

#### Neural network pruning.

As neural networks turn to become deeper and more sophisticated, model pruning
technology has drawn great attention than ever in the last decade, as less
parameters will facilitate deep model deployment in practical applications
[han2015deep, denil2013overparamterized_and_redundancy]. With the goal of
finding high-sparse, high-accuracy subnetworks from original dense models, a
variety of pruning methods have been developed such as heuristics-based
pruning [zhang2021efficient, frankle2018lottery, wang2020pick, lee2018snip,
tanaka2020pruning, alizadeh2022prospect] and optimization-based pruning
[sehwag2020hydra, lee2020layer, ramanujan2020sHiddenSubnetwork,
csordas2020neural]. The former identifies redundant model weights by
leveraging heuristics-based metrics such as magnitudes [janowsky1989pruning,
han2015learning] SL: Add more from Jonathan Frankle’s work and Tianlong’s work
, gradient and Hessian statistics [molchanov2019importance, yao2020pyhessian]
SL: [add grasp, SynFlow, snip?? do they belong to?], Taylor coefficients
[mozer1989skeletonization, lecun1990optimal, hassibi1993second,
molchanov2016pruning] SL: [I did not understand Taylor coefficients. Please be
more specific and preciser.]. The latter is typically built on sparsity-
promoting optimization [liu2017learning, he2017channel, zhou2016less,
ren2018admmnn, guo2021gdp, sehwag2020hydra] SL: [add admm pruning from
Yanzhi’s group; no need to cite optimization literature], where model weights
are training by penalizing their sparsity-inducing norms, such as $\ell_{1}$
norm for irregular weight pruning and $\ell_{2}$ norm for structured pruning.

In addition, based on whether or not the targeted model sparsity is achieved
using a one-shot or iterative pruning schedule, existing methods can be
further split into one-shot pruning and iterative pruning. Among them, the
Iterative Magnitude Pruning (IMP) [frankle2018lottery, gale2019state,
pmlr-v139-zhang21c, chen2020lottery, chen2020lottery2, yu2019playing,
chen2021gans, ma2021good, gan2021playing, chen2021unified,
kalibhat2021winning, chen2021gans, chen2021ultra] has played a significant
role in identifying high-quality sub-networks, as supported by the Lottery
Ticket Hypothesis (LTH) [Renda2020Comparing, frankle2020linear], which unveils
the existence of ‘winning tickets’, i.e., sub-networks that can be trained in
isolation to reach or even succeed the performance of the original dense
model. However, it is computationally expensive to find winning tickets by IMP
due to its iterative nature. As a computation-efficient counterpart of IMP,
One-shot Magnitude Pruning (OMP) [frankle2018lottery] is developed, where a
dense model is directly pruned to the targeted sparsity. The other one-shot
pruning methods such as pruning at initialization [wang2020pick, lee2018snip,
tanaka2020pruning, alizadeh2022prospect] further reduce the computation cost
of pruning by avoiding the dense model training. Yet, it is worth noting that
the gain of efficiency using one-shot pruning is accompanied by the loss of
pruning accuracy, compared to iterative pruning.

Most of the aforementioned literature focuses on irregular weight pruning,
which regards the weight sparsity number as the main pruning purpose. By
contrast, structured pruning takes into consideration the targeted sparse
patterns, such as filter, channel, and kernel-wise pruning [he2017channel,
li2016pruning, guo2021gdp, chen2022coarsening]. SL: [please cite more e.g.,
pattern-aware pruning from yanzhi’s group]. The structured pruning method is
especially favored for actual deep model deployment at the presence of
hardware implementation constraints SL: [refs]. However, compared to irregular
weight pruning, it imposes challenges to preserve pruning accuracy and thus
becomes difficult to obtain structure-aware winning tickets SL: [refs].

SL: [I stop here. @Pari, if possible, please feel free to check the following.
I will do it tonight as well.] To bridge the gap between the one-shot pruning
and IMP and explore the possibility of achieving winning tickets in structured
pruning, we revisit pruning problem from the perspective of bi-level
optimization. In our work, we for the first time propose an efficient first-
order pruning framework through bi-level optimization for both structured and
unstructured pruning. In particular, our proposed algorithm solves the bi-
level pruning algorithm and achieves better performance in various settings
than IMP with computational costs comparable to OMP.

#### Bi-level optimization (BLO).

BLO is a general hierarchical optimization framework, where the upper-level
problem highly dependent on the solution to the lower-level problem. Such a
nested structure makes the BLO in its most generic form very difficult to
solve, since it is hard to characterize the small influence from the lower-
level optimizer on the upper-level problem. Thus, existing works devoted to
the design of algorithms and theory for various special cases of BLO problems.
Applications range from the classical approximate descent methods
[falk1995bilevel, vicente1994descent], penalty-based methods
[white1993penalty, gould2016differentiating], to recent BigSAM
[sabach2017first] and its extensions [liu2020generic, li2020improved]. It is
also studied in the area of stochastic algorithms [ghadimi2018approximation,
ji2020bilevel, hong2020two] and back-propagation based algorithms
[franceschi2017forward, grazzi2020iteration, shaban2019truncated]. Besides,
BLO has also advanced adversarial training [zhang2021bat], meta-learning
[rajeswaran2019meta], data poisoning attack generation [huang2020metapoison],
network architecture search [liu2018darts] as well as reinforcement learning
[chen2019research]. However, to the best of our knowledge, none of the
existing works adopt the BLO approach for neural network pruning tasks. Our
work, for the first time, advance the algorithmic foundation of neural network
pruning through the lens of BLO.

##  3 Preliminaries and Motivation

In this section, we begin by providing a brief background of neural network
pruning, including the predominant iterative magnitude pruning (IMP) method
and the associated lottery ticket hypothesis (LTH) [frankle2018lottery] to
find high-quality subnetworks (known as ‘winning tickets’) that–when trained
in isolation–reach test accuracy comparable to or even better than the
original dense network. However, the superior performance of IMP is
accompanied with a prohibitive computation cost, as will be exemplified in
Fig. 2. In order to overcome the computational difficulty of IMP, we next
introduce the one-shot pruning method and elaborate on its comparison with
IMP. Spurred by the pros and cons of current pruning methods, we eventually
state the problem to be studied in this paper.

#### Heuristics-based IMP.

The heuristics-based IMP is the currently prevailing approach to identify the
‘winning ticket’ in model pruning, so as to attain model sparsity without
suffering performance loss [frankle2018lottery, Renda2020Comparing,
frankle2020linear, blalock2020state, chen2020lottery, chen2021lottery]. For
example, the current state-of-the-art (SOTA) implementation of IMP in
[ma2021sanity] can lead to a pruned ResNet-20 on CIFAR-10 with 74% sparsity
and $92.12\%$ test accuracy, which is comparable to the performance of the
original dense model; see Fig. 2. The algorithm of IMP typically contains two
key ingredients: (i) a temporally-evolved pruning schedule (given by the
growing model sparsity over pruning iterations), and (ii) the pruning-
retraining learning mechanism applied at each pruning round. For (i),
supposing that there exist $n$ pruning rounds, each round prunes
$(p\%)^{\frac{1}{n}}$ of the model weights that survive the previous round,
where $p\%$ denotes the final pruning ratio to be achieved. For (ii), the non-
zero remaining weights at each pruning round are re-set to the weights at
initialization or at an early-training epoch [Renda2020Comparing], and re-
trained till convergence. In brief, IMP repeatedly trains, prunes, and resets
the network over multiple rounds.

Figure 2:  Illustration of the pros and cons of different pruning methods
executed over (ResNet-20, CIFAR-10). (a) Test accuracy vs. different pruning
ratios: IMP significantly outperforms one-shot pruning methods (OMP
[frankle2018lottery], ProsPr [alizadeh2022prospect], Grasp [wang2020pick],
Snip [lee2018snip], SynFlow [tanaka2020pruning]). (b) Computation time of
pruning to different sparsity levels: IMP becomes much more computational
intensive when pruning is targeted for a higher sparsity level.

Yet, the superior pruning accuracy of IMP is at the cost of computation
efficiency. For ease of illustration, Fig. 2 shows the example of pruning
ResNet-20 over CIFAR-10 and demonstrates the accuracy and the computation time
of various pruning methods versus the growing sparsity ratio. As we can see,
IMP outperforms the other pruning baselines and finds the sparsest winning
ticket (marked by the red star in Fig. 2a) using $n=6$ pruning rounds.
However, this leads to $3\times$ more than the cost of training a dense model
from scratch.

#### IMP vs. one-shot pruning.

To reduce the computation cost of IMP, different kinds of ‘accelerated’
pruning methods were developed [zhang2021efficient, frankle2018lottery,
wang2020pick, lee2018snip, tanaka2020pruning, alizadeh2022prospect]; And many
of them fall into the one-shot pruning category: The network is directly
pruned to the eventual target sparsity level and retrained once. In
particular, the one-shot magnitude pruning (OMP) has been served as an
important baseline that simplifies IMP directly [frankle2018lottery]. OMP
obeys the pruning-retraining paradigm, but it prunes the model to the target
sparsity using only one pruning round. However, the promise of current one-
shot pruning methods is dataset/model-specific [ma2021sanity,
su2020sanity_random_win] and mostly lie in the low-sparsity regime
[alizadeh2022prospect]. As systematically studied in [frankle2020pruning,
ma2021sanity], there exists a clear performance gap between one-shot pruning
and IMP. Fig. 2 also justifies the above limitation. As we can see, although
one-shot pruning is computationally light (Fig. 2b), it suffers from an
evident accuracy drop compared to IMP (Fig. 2a), especially at the large
sparsity levels. Thus, there exists a strong demand for the algorithmic
advancement of model pruning that can enjoy the dual advantages of one-shot
pruning and IMP.

#### Problem Statement.

As described above, the existing one-shot pruning methods do not suffice to
close the performance gap with IMP. One possible root cause is its strong
dependence on the heuristics-based weight sparsification strategy. This leaves
the question of whether computation efficiency and pruning accuracy could be
co-improved through a more rigorous and dedicated optimization-based approach.
Thus, we ask: Beyond heuristics, how to establish the optimization basis of a
successful pruning mechanism that can attain high pruning accuracy (like IMP)
as well as high computation efficiency (like one-shot pruning)?

##  4 BiP: Model Pruning via Bi-level Optimization

In this section, we re-investigate the problem of model pruning through the
lens of bi-level optimization (BLO), and show that BLO not only provides a
natural interpretation of the pruning-retraining learning paradigm but also
brings in additional design flexibilities for model pruning. Next, we
establish the algorithmic foundation of BLO-enabled model pruning (termed bi-
level pruning or BiP), extended from the implicit gradient approach. Further,
we show that BiP can be solved as easily as the first-order alternating
optimization by taking advantages of the bi-linear nature of the pruning
variables. Lastly, we summarize the BiP algorithm and its implementation
details.

#### A BLO viewpoint on model pruning

As described in Sec. 3, the pruning-retraining learning paradigm covers two
kinds of tasks: ❶ pruning to determine the sparse pattern of model weights,
and ❷ training remaining non-zero weights to recover the model accuracy. In
existing optimization-based pruning methods [wen2016learning,
zhang2018systematic, hoefler2021sparsity, wang2020proximal], the tasks ❶-❷ are
achieved by merely optimizing model weights, together with penalizing their
sparsity-inducing norms, e.g., the $\ell_{1}$ and $\ell_{2}$ norms
[bach2012optimization]. Different from the above formulation, we propose to
separate optimization variables involved in the pruning tasks ❶ and ❷. This
leads to the (binary) pruning mask variable $\mathbf{m}\in\\{0,1\\}^{n}$ and
the model weight variable $\boldsymbol{\theta}\in\mathbb{R}^{n}$,
respectively. Here $n$ denotes the total number of model parameters.
Accordingly, the pruned model is given by
$(\mathbf{m}\odot\boldsymbol{\theta})$, where $\odot$ denotes the element-wise
multiplication. As will be evident later, this variable disentanglement
enables us to explicitly depict the co-evolving of pruning and retraining
through BLO, and helps us to customize the pruning task with high flexibility.

In what follows, we elaborate on how BLO can be established to co-optimize the
pruning mask $\mathbf{m}$ and the retrained sparse model weights
$\boldsymbol{\theta}$. Given the pruning ratio $p\%$, the sparsity constraint
is given by $\mathbf{m}\in\mathcal{S}$, where
$\mathcal{S}=\\{\mathbf{m}\,|\,\mathbf{m}\in\\{0,1\\}^{n},\mathbf{1}^{T}\mathbf{m}\leq
k\\}$ and $k=p\%n$. Our goal is to prune the original model directly to the
targeted sparse number $k$ (i.e., without calling for the IMP-like sparsity
schedule as described in Sec. 3) and obtain the optimized sparse model
$(\mathbf{m}\odot\boldsymbol{\theta})$. To this end, we interpret the pruning
task (i.e., ❶) and the model retraining task (i.e., ❷) as two optimization
levels, where the former is formulated as an upper-level optimization problem,
and it relies on the optimizer of the lower-level retraining task. We then
cast the model pruning problem as the following BLO problem (with ❷ being
nested inside ❶): 55todo: 5Mention structured pruning? For structured pruning,
manipulate m with Kronecker and transfer the element-wise pruning unit to a
structured one.

$\displaystyle\displaystyle\begin{array}[]{llll}\displaystyle\operatorname*{\text{minimize}}_{\mathbf{\mathbf{m}}\in\mathcal{S}}&\underbrace{\ell(\mathbf{m}\odot{\boldsymbol{\theta}}^{*}(\mathbf{m}))}_{\text{
\char 182: Pruning task}};&\quad\operatorname*{\text{subject
to}}&\underbrace{\displaystyle{\boldsymbol{\theta}}^{*}(\mathbf{m})=\argmin_{\boldsymbol{\theta}\in\mathbb{R}^{n}}\overbrace{\ell(\mathbf{m}\odot\boldsymbol{\theta})+\frac{\gamma}{2}\|\boldsymbol{\theta}\|_{2}^{2}}^{\text{$\mathrel{\mathop{:}}=g(\mathbf{m},\boldsymbol{\theta})$}}}_{\text{\char
183: Sparsity-fixed model retraining}},\end{array}$ (2)

where $\ell$ denotes the model training loss (e.g., cross entropy),
$\mathbf{m}$ and $\boldsymbol{\theta}$ are the upper-level and lower-level
optimization variables respectively, ${\boldsymbol{\theta}}^{*}(\mathbf{m})$
signifies the lower-level solution by minimizing the objective function $g$
given the pruning mask $\mathbf{m}$, and $\gamma>0$ is a regularization
parameter introduced to strongly convexify the lower-level optimization so as
to stabilize the gradient flow from $\boldsymbol{\theta}^{*}(\mathbf{m})$ to
$\mathbf{m}$ and thus the convergence of BLO [shaham2015understanding,
hong2020two]. In a sharp contrast to existing single-level optimization-based
model pruning methods [wen2016learning, zhang2018systematic,
hoefler2021sparsity, wang2020proximal], the BLO formulation (2) brings in two
advantages.

First, BLO has flexibility to use mismatched pruning and retraining objectives
at the upper and lower optimization levels, respectively. This flexibility
exists for customizing not only the regularized lower-level objective function
in (2) but also optimization algorithm per level. To be more specific, one can
update the upper-level pruning mask $\mathbf{m}$ using a data batch (called
$\mathcal{B}_{2}$) distinct from the one (called $\mathcal{B}_{1}$) used for
obtaining the lower-level solution ${\boldsymbol{\theta}}^{*}(\mathbf{m})$.
The resulting BLO procedure can then mimic the idea of meta learning to
improve model generalization [deng2021meta], where the lower-level problem
fine-tunes $\boldsymbol{\theta}$ using $\mathcal{B}_{1}$, and the upper-level
problem validates the generalization of the sparsity-aware finetuned model
$(\mathbf{m}\odot\boldsymbol{\theta}^{*}(\mathbf{m}))$ using
$\mathcal{B}_{2}$.

Second, BLO enables us to explicitly model and optimize the coupling between
the retained model weights $\boldsymbol{\theta}^{*}(\mathbf{m})$ and the
pruning mask $\mathbf{m}$ through the implicit gradient (IG)-based
optimization routine. Here IG refers to the gradient of the lower-level
solution $\boldsymbol{\theta}^{*}(\mathbf{m})$ with respect to (w.r.t.) the
upper-level variable $\mathbf{m}$, and its derivation calls the implicit
function theory [gould2016differentiating]. The use of IG makes our proposed
BLO-oriented pruning (2) significantly different from the greedy alternating
minimization [csiszar1984information] that learns the upper-level and lower-
level variables independently (i.e., minimizes one variable by fixing the
other). We refer readers to the following section for the detailed IG theory.
We will also show in Sec. 5 that the pruning strategy acquired from (2) can
outperform IMP in many pruning scenarios but is much more efficient in
computation as it does not call for the scheduler of iterative pruning ratios.

#### Optimization foundation of BiP.

The key optimization challenge of solving the BiP problem (2) lies in the
computation of IG (implicit gradient). Prior to developing an effective
solution, we first elaborate on the IG challenge, the unique characteristic of
BLO. In the context of gradient descent, the gradient of the objective
function in (2) yields

$\displaystyle\displaystyle\underbrace{\frac{d\ell(\mathbf{m}\odot\boldsymbol{\theta}^{*}(\mathbf{m}))}{d\mathbf{m}}}_{\text{Gradient
of
objective}}=\nabla_{\mathbf{m}}\ell(\mathbf{m}\odot{\boldsymbol{\theta}}^{*}(\mathbf{m}))+\underbrace{\frac{{d({\boldsymbol{\theta}}^{*}(\mathbf{m})^{\top})}}{d\mathbf{m}}}_{\text{{IG}}}\nabla_{\boldsymbol{\theta}}\ell(\mathbf{m}\odot\boldsymbol{\theta}^{*}(\mathbf{m})),$
(3)

where $\nabla_{\mathbf{m}}$ and $\nabla_{\boldsymbol{\theta}}$ denote the
partial derivatives of the bi-variate function
$\ell(\mathbf{m},\boldsymbol{\theta})$ w.r.t. the variable $\mathbf{m}$ and
$\boldsymbol{\theta}$ respectively,
$d\boldsymbol{\theta}^{\top}/d\mathbf{m}\in\mathbb{R}^{n\times n}$ denotes the
vector-wise full derivative, and for ease of notation, we will omit the
transpose $\top$ when the context is clear. In (3), the IG challenge refers to
the demand for computing gradients of the implicit function
$\boldsymbol{\theta}^{*}(\mathbf{m})=\argmin_{\boldsymbol{\theta}}g(\mathbf{m},\boldsymbol{\theta})$
w.r.t. $\mathbf{m}$, where recall from (2) that
$g(\mathbf{m},\boldsymbol{\theta})\mathrel{\mathop{:}}=\ell(\mathbf{m}\odot\boldsymbol{\theta})+\frac{\gamma}{2}\|\boldsymbol{\theta}\|_{2}^{2}$.

Next, we derive the IG formula following the more rigorous implicit function
theory [gould2016differentiating, hong2020two, rajeswaran2019meta]. Based on
the fact that the $\boldsymbol{\theta}^{*}(\mathbf{m})$ is satisfied with the
stationary condition of the lower-level objective function in (3), it is not
difficult to obtain that

$\displaystyle\displaystyle\frac{d\boldsymbol{\theta}^{*}(\mathbf{m})}{d\mathbf{m}}=-\nabla^{2}_{\mathbf{m}\boldsymbol{\theta}}\ell(\mathbf{m}\odot\boldsymbol{\theta}^{*})[\nabla^{2}_{\boldsymbol{\theta}}\ell(\mathbf{m}\odot\boldsymbol{\theta}^{*})+\gamma\mathbf{I}]^{-1},$
(4)

where $\nabla_{\mathbf{m}\boldsymbol{\theta}}^{2}\ell$ and
$\nabla_{\boldsymbol{\theta}\boldsymbol{\theta}}^{2}\ell$ denote the second-
order partial derivatives of $\ell$ respectively, and $(\cdot)^{-1}$ denotes
the matrix inversion operation. We refer readers to Appendix LABEL:eq:_app_IG
for the full derivation of (4).

Yet, the exact IG formula (4) remains difficult in calculation due to the
presence of matrix inversion and second-order partial derivatives. 66todo:
6$\gamma$ is critical in this derivation. Without $\gamma$, and with the
$\nabla^{2}_{\boldsymbol{\theta}}=0$ assumption, this does not go through. How
do we best motivate the introduction of $\gamma$ in (1) since standard pruning
setups and schemes dont introduce such regularization term? To simplify it, we
impose the Hessian-free assumption,
$\nabla^{2}_{\boldsymbol{\theta}}\ell=\mathbf{0}$, which is mild in general;
For example, the decision boundaries of neural networks with ReLU activations
are piece-wise linear in a tropical hyper-surface [alfarra2020decision], and
this assumption has been widely used in meta learning [finn2017model]. Given
$\nabla^{2}_{\boldsymbol{\theta}}\ell=\mathbf{0}$, the matrix inversion in (4)
can be then mitigated, leading to the IG formula

$\displaystyle\frac{d\boldsymbol{\theta}^{*}(\mathbf{m})}{d\mathbf{m}}=-\frac{1}{\gamma}\nabla^{2}_{\mathbf{m}\boldsymbol{\theta}}\ell(\mathbf{m}\odot\boldsymbol{\theta}^{*}).$
(5)

At the first glance, the computation of the simplified IG (5) still requires
the mixed (second-order) partial derivative
$\nabla^{2}_{\mathbf{m}\boldsymbol{\theta}}\ell$. However, BiP is a special
class of BLO problems with bi-linear variables
$(\mathbf{m}\odot\boldsymbol{\theta})$. Based on this bi-linearity, we can
prove that IG in (5) can be analytically expressed using only first-order
derivatives; see the following theorem.

###### Proposition 1

Assuming $\nabla^{2}_{\boldsymbol{\theta}}\ell=0$, the implicit gradient (5)
is then given by

$\displaystyle\displaystyle\frac{d\boldsymbol{\theta}^{*}(\mathbf{m})}{d\mathbf{m}}=-\frac{1}{\gamma}\mathrm{diag}(\nabla_{\mathbf{z}}\ell(\mathbf{z})\left|{}_{\mathbf{z}={\mathbf{m}}\odot{\boldsymbol{\theta}}^{*}}\right.);$
(6)

Further, the gradient of the objective function given by (3) becomes

$\displaystyle\displaystyle\frac{d\ell(\mathbf{m}\odot\boldsymbol{\theta}^{*})}{d\mathbf{m}}=(\boldsymbol{\theta}^{*}-\frac{1}{\gamma}\mathbf{m}\odot\nabla_{\mathbf{z}}\ell(\mathbf{z}))\odot\nabla_{\mathbf{z}}\ell(\mathbf{z}),$
(7)

where for ease of notation $\nabla_{\mathbf{z}}\ell(\mathbf{z})$ signifies
$\nabla_{\mathbf{z}}\ell(\mathbf{z})\left|{}_{\mathbf{z}={\mathbf{m}}\odot{\boldsymbol{\theta}}^{*}}\right.$,
and recall that $\odot$ denotes the element-wise multiplication.

Proof: Using chain-rule, we can obtain that

$\displaystyle\nabla_{\boldsymbol{\theta}}\ell(\mathbf{m}\odot\boldsymbol{\theta}^{*})=\mathrm{diag}(\mathbf{m})\nabla_{\mathbf{z}}\ell(\mathbf{z})=\mathbf{m}\odot\nabla_{\mathbf{z}}\ell(\mathbf{z});$
(8) similarly,
$\displaystyle\nabla_{\mathbf{m}}\ell(\mathbf{m}\odot\boldsymbol{\theta}^{*})=\mathrm{diag}(\boldsymbol{\theta}^{*})\nabla_{\mathbf{z}}\ell(\mathbf{z})=\boldsymbol{\theta}^{*}\odot\nabla_{\mathbf{z}}\ell(\mathbf{z})$
(9)

where $\mathrm{diag}(\cdot)$ represents a diagonal matrix with $\cdot$ being
the main diagonal vector. Further, we can convert (5) to

$\displaystyle\displaystyle\nabla^{2}_{{\mathbf{m}}{\boldsymbol{\theta}}}\ell({\mathbf{m}}\odot{\boldsymbol{\theta}^{*}})$
$\displaystyle\overset{\eqref{eq:
chainrule_theta}}{=}\nabla_{{\mathbf{m}}}\left[{\mathbf{m}}\odot\nabla_{\mathbf{z}}\ell(\mathbf{z})\right]\overset{\text{chain
rule}}{=}\mathrm{diag}({\nabla_{\mathbf{z}}\ell(\mathbf{z})})+\mathrm{diag}({\mathbf{m}})[\nabla_{{\mathbf{m}}}(\nabla_{\mathbf{z}}\ell(\mathbf{z}))]$
$\displaystyle\overset{\eqref{eq:
chainrule_mask}}{=}\mathrm{diag}({\nabla_{\mathbf{z}}\ell(\mathbf{z})})+\mathrm{diag}({\mathbf{m}})[\mathrm{diag}(\boldsymbol{\theta}^{*})\nabla^{2}_{\mathbf{z}}\ell(\mathbf{z})]=\mathrm{diag}(\nabla_{\mathbf{z}}\ell(\mathbf{z})),$
(10)

where the last equality holds due to the Hessian-free assumption. With (10)
and (5) we can prove (6).

Next, substituting the IG (6) to the upper-level gradient (3), we obtain that

$\displaystyle\displaystyle\frac{d\ell(\mathbf{m}\odot\boldsymbol{\theta}^{*})}{d\mathbf{m}}$
$\displaystyle\,\,\,\,=\nabla_{\mathbf{m}}\ell({\mathbf{m}}\odot{\boldsymbol{\theta}}^{*})-{\frac{1}{\gamma}\nabla_{\mathbf{z}}\ell(\mathbf{z})\odot\nabla_{{\boldsymbol{\theta}}}\ell({\mathbf{m}}\odot{\boldsymbol{\theta}^{*}})}$
$\displaystyle\underset{}{\overset{\eqref{eq: chainrule_theta},\eqref{eq:
chainrule_mask}}{=}}\boldsymbol{\theta}^{*}\odot\nabla_{\mathbf{z}}\ell(\mathbf{z})-\frac{1}{\gamma}\nabla_{\mathbf{z}}\ell(\mathbf{z})\odot(\mathbf{m}\odot\nabla_{\mathbf{z}}\ell(\mathbf{z}))=(\boldsymbol{\theta}^{*}-\frac{1}{\gamma}\mathbf{m}\odot\nabla_{\mathbf{z}}\ell(\mathbf{z}))\odot\nabla_{\mathbf{z}}\ell(\mathbf{z}),$

which leads to (7). The proof is now complete. $\square$ 77todo: 7Do we have
any empirical evidence that shows that using (7) is better than using (9)?
That would be good to highlight advantage of BLO and IG. YY: Pari means the
significance of the extra term brought by IG.

The key insight drawn from Prop. 1 is that the bi-linearity of pruning
variables (i.e., $\mathbf{m}\odot\boldsymbol{\theta}^{*}$) makes the IG-
involved gradient (3) easily computed, and the computational complexity is
almost same as that of computing the first-order gradient
$\nabla_{\mathbf{z}}\ell(\mathbf{z})$ just once, as supported by (7).

###### Corollary 1

Given $\mathbf{m}\in\\{0,1\\}^{n},\boldsymbol{\theta}\in\mathbb{R}^{n}$ and
that the computation of the gradient
$\nabla_{\mathbf{z}}\ell(\mathbf{z})\in\mathbb{R}^{n}$ takes $O(n)$ time, the
computation of the loss gradient (3) with respect to the mask $\mathbf{m}$ of
the objective (2) using (7) in Proposition 1 takes $O(n)$ time.

####  BiP algorithm and implementation.

We will next formalize the BiP algorithm based on Prop. 1 and the alternating
gradient descent method [hong2020two]. At iteration $t$, there exist two main
steps.

88todo: 8YY: Hi Pari, feel free to elaborate on the algorithm with the way
that you think is the best

★ Lower-level SGD for model retraining: Given $\mathbf{m}^{(t-1)}$,
$\boldsymbol{\theta}^{(t-1)}$, and
$\mathbf{z}^{(t-1)}\mathrel{\mathop{:}}=\mathbf{m}^{(t-1)}\odot\boldsymbol{\theta}^{(t-1)}$,
we update $\boldsymbol{\theta}^{(t)}$ by applying SGD (stochastic gradient
descent) to the lower-level problem of (2),

$\displaystyle\boldsymbol{\theta}^{(t)}$
$\displaystyle=\boldsymbol{\theta}^{(t-1)}-\alpha\nabla_{\boldsymbol{\theta}}g(\mathbf{m}^{(t-1)},\boldsymbol{\theta}^{(t-1)})$
$\displaystyle\overset{\eqref{eq:
chainrule_theta}}{=}\boldsymbol{\theta}^{(t-1)}-\alpha[\mathbf{m}^{(t-1)}\odot\nabla_{\mathbf{z}}\ell(\mathbf{z})\left|\right._{\mathbf{z}=\mathbf{z}^{t-1}}+\gamma\boldsymbol{\theta}^{(t-1)}],$
($\boldsymbol{\theta}$-step)

where $\alpha>0$ is the lower-level learning rate.

★ Upper-level PGD for pruning: Given $\mathbf{m}^{(t-1)}$,
$\boldsymbol{\theta}^{(t)}$, and
$\mathbf{z}^{(t+1/2)}\mathrel{\mathop{:}}=\mathbf{m}^{(t-1)}\odot\boldsymbol{\theta}^{(t)}$,
we update $\mathbf{m}$ using PGD (projected gradient descent) along the IG-
enhanced descent direction (3),

$\displaystyle\mathbf{m}^{(t)}=$
$\displaystyle\mathcal{P}_{\mathcal{S}}\left[\mathbf{m}^{(t-1)}-\beta\frac{d\ell(\mathbf{m}\odot\boldsymbol{\theta}^{(t)})}{d\mathbf{m}}\left|\right._{\mathbf{m}=\mathbf{m}^{(t-1)}}\right]$
$\displaystyle\overset{\eqref{eq: upper_update_final}}{=}$
$\displaystyle\mathcal{P}_{\mathcal{S}}\left[\mathbf{m}^{(t-1)}-\beta\left(\boldsymbol{\theta}^{(t)}-\frac{1}{\gamma}\mathbf{m}^{(t-1)}\odot\nabla_{\mathbf{z}}\ell(\mathbf{z})\left|\right._{\mathbf{z}=\mathbf{z}^{(t+1/2)}}\right)\odot\nabla_{\mathbf{z}}\ell(\mathbf{z})\left|\right._{\mathbf{z}=\mathbf{z}^{(t+1/2)}}\right],$
($\mathbf{m}$-step)

where $\beta>0$ is the upper-level learning rate, and
$\mathcal{P}_{\mathcal{S}}(\cdot)$ denotes the Euclidean projection onto the
constraint set $\mathcal{S}$ given by
$\mathcal{S}=\\{\mathbf{m}\,|\,\mathbf{m}\in\\{0,1\\}^{n},\mathbf{1}^{T}\mathbf{m}\leq
k\\}$ in (2) and is achieved by the top-$k$ hard-thresholding operation as
will be detailed later.

In BiP, the ($\boldsymbol{\theta}$-step) and ($\mathbf{m}$-step) steps execute
iteratively till convergence. In practice, we find that BiP converges very
well; see Fig. LABEL:fig:_convergence for an ablation study. It is also worth
noting that both ($\boldsymbol{\theta}$-step) and ($\mathbf{m}$-step) steps
only require the first-order information
$\nabla_{\mathbf{z}}\ell(\mathbf{z})$, demonstrating that the BiP algorithm
can be conducted as efficiently as first-order optimization. We present more
implementation details of BiP below.

✦ Discrete optimization over $\mathbf{m}$: We follow the ‘convex relaxation +
hard thresholding’ mechanism used in [sehwag2020hydra,
ramanujan2020sHiddenSubnetwork]. Specifically, we relax the binary masking
variables to continuous masking scores $\mathbf{m}\in[\mathbf{0},\mathbf{1}]$.
We then acquire loss gradients at the backward pass based on the relaxed
$\mathbf{m}$. At the forward pass, we project it onto the discrete constraint
set $\mathcal{S}$ using the hard thresholding operator, where its top $k$
elements are set by $1$s and the others by $0$s.

✦ Data batch selection for lower-level and upper-level optimization: We adopt
different data batches (with the same batch size) when implementing
($\boldsymbol{\theta}$-step) and ($\mathbf{m}$-step). This is one BLO’s
advantage of, which owns the flexibility to customize the lower-level and
upper-level optimization problems. And the use of diverse data batches is
beneficial to model generalization as shown in [deng2021meta].

✦ Hyperparameter tuning: As described in
($\boldsymbol{\theta}$-step)-($\mathbf{m}$-step), BiP needs to set two
learning rates $\alpha$ and $\beta$ for lower-level and upper-level
optimization, respectively. We choose $\alpha=0.01$ and $\beta=0.1$ in all
experiments, where we adopt the mask learning rate $\beta$ from Hydra
[sehwag2020hydra] and set a smaller lower-level learning rate $\alpha$, as
$\boldsymbol{\theta}$ is initialized by a pre-trained dense model. We show
ablation study on $\alpha$ in Fig. 7(c). BLO also brings in the low-level
convexification parameter $\gamma$. We set $\gamma=1.0$ in experiments and
refer readers to Fig. 7(b) for a sanity check.

✦ One-shot vs. multi-step SGD: In ($\boldsymbol{\theta}$-step), the one-step
SGD is used and helps reduce the computation overhead. In practice, we also
find that the one-step SGD is sufficient: The use of multi-step SGD in BiP
does not yield much significant improvement over the one-step version; see
Fig. 7(a).

✦ Extension to structured pruning: We formulate and solve the BiP problem in
the context of unstructured (element-wise) weight pruning. However, if define
the pruning mask $\mathbf{m}$ corresponding to model’s structure units (e.g.,
filters), then it is easy to apply BiP to structured weight pruning; see Fig.
4 and Fig. LABEL:fig:_structured_performance_cifar_appendix for empirical
results.

##  5 Experiments

In this section, we present extensive experiment results to show the
effectiveness of BiP across multiple model architectures, various datasets,
and different pruning setups. Compared to IMP, one-shot pruning, and
optimization-based pruning baselines, we find that BiP leads to better winning
tickets in most cases and is computationally as efficient as one-shot pruning.

###  5.1 Experiment Setup

#### Datasets and models.

Following the pruning benchmark in [ma2021sanity], we consider 4 datasets
including CIFAR-10 [krizhevsky2009learning], CIFAR-100
[krizhevsky2009learning], Tiny-ImageNet [le2015tiny], ImageNet
[deng2009imagenet], and 5 architecture types including ResNet-20/56/18/50 and
VGG-16 [he2016deep, simonyan2014very]. Tab. LABEL:table:_impl_details
summarizes these dataset and model information and setups.

#### Baselines, training, and evaluation.

When comparing with BiP, we mainly focus on 4 state-of-the-art (SOTA) pruning
methods, ① IMP [frankle2018lottery], ② OMP [frankle2018lottery], ③ Grasp
[wang2019picking] (a one-shot pruning method by analyzing gradient flow at
initialization), and ④ Hydra [sehwag2020hydra] (an optimization-based pruning
method that optimizes masking scores). It is worth noting that there exist
various implementations of IMP, e.g., specified by different learning rate and
model initialization ‘rewinding’ strategies [Renda2020Comparing]. To make a
fair comparison, we follow the recent IMP benchmark in [ma2021sanity], which
can find the best winning tickets over current heuristics-based pruning
baselines. We also remark that Hydra is originally proposed for improving
adversarial robustness of a pruned model, but it can be easily customized for
standard pruning when setting the adversary’s strength as $0$
[sehwag2020hydra]. We choose Hydra as a baseline because it can be regarded as
a specification of BiP, given by a single-level sparse mask optimization task
together with post-optimization weight re-training. When implementing BiP,
unless specified otherwise, we use the $1$-step SGD in
($\boldsymbol{\theta}$-step), and set the learning rates $(\alpha,\beta)$ and
the lower-level regularization parameter $\gamma$ as described in the previous
section. When implementing baselines, we follow their official code repository
setups. We evaluate the performance of BiP and the above baselines mainly from
two perspectives:14 (1) the test accuracy of the obtained sub-network, and (2)
the runtime of pruning to reach the desired sparsity. We refer readers to Tab.
LABEL:tab:_training_details, Appendix LABEL:app:_training for more training
and evaluation details, such as training epochs, learning rate schedulers, and
so forth.

###  5.2 Experiment Results

####  BiP identifies high-accuracy winning tickets.

In what follows, we peer into the quality of winning tickets identified by
BiP. Two key observations can be drawn from our experiment results: (1) BiP
finds winning tickets of higher accuracy and/or higher sparsity than the
baselines in most cases (as shown in Fig. 3 and Tab. 1; (2) The superiority of
BiP holds for both unstructured pruning and structured pruning (as shown in
Fig. 4 and Fig. LABEL:fig:_structured_performance_cifar_appendix ).

Figure 3: Unstructured pruning trajectory given by test accuracy (%) vs.
sparsity (%) on various dataset-model pairs. The proposed BiP is compared with
Hydra [sehwag2020hydra], IMP [ma2021sanity], OMP [ma2021sanity], Grasp
[wang2019picking]. And the performance of dense model and that of the best
winning ticket are marked using dashed lines in each plot. The solid line and
shaded area of each pruning method represent the mean and variance of test
accuracies over 3 independent trials. We observe that BiP consistently
outperforms the other baselines.  Table 1: The sparsest winning tickets found
by different methods at different data-model setups. YY: Winning tickets refer
to the sparse models with average test accuracy no less than the dense ones.
In each cell, $p\%$ (acc$\pm$std%) represents the sparsity as well as the test
accuracy. The test accuracy of dense models can be found in the header. ✗
signifies that no winning ticket is found by a pruning method. Given the data-
model setup (i.e., per column), the sparsest winning ticket is highlighted in
bold. YY: test accuracy and std

Method | CIFAR-10 | CIFAR-100
---|---|---
ResNet-18 | ResNet-20 | ResNet-56 | VGG-16 | ResNet-18 | ResNet-20
(94.77%) | (92.12%) | (92.95%) | (93.65%) | (76.67%) | (68.69%)
IMP | 87% (94.77$\pm$0.10%) |  74% (92.15$\pm$0.15%) |  74% (92.99$\pm$0.12%) | 89% (93.68$\pm$0.05%) | 87% (76.91$\pm$0.19%) | ✗
OMP | 49% (94.80$\pm$0.10%) | ✗ | ✗ | 20% (93.79$\pm$0.06%) | 74% (76.99$\pm$0.07%) | ✗
Grasp | ✗ | ✗ | 36% (93.07$\pm$0.34%) | ✗ | ✗ | ✗
Hydra | ✗ | ✗ | ✗ | 87% (93.73$\pm$0.03%) | ✗ | 20% (68.94$\pm$0.17%)
BiP |  89% (94.79$\pm$0.15%) | 67% (92.14$\pm$0.15%) |  74% (93.13$\pm$0.04%) |  93% (93.75$\pm$0.15%) |  89% (76.69$\pm$0.18%) |  49% (68.78$\pm$0.10%)

Figure 4: Filter pruning trajectory given by test accuracy (%) vs. pruning
ratio (%). The visual presentation setting is the same as Fig. 3. We observe
that BiP identifies winning tickets of structured pruning in certain sparsity
regimes.

Fig. 3 shows the unstructured pruning trajectory (given by test accuracy vs.
sparsity ratio) of BiP and baseline methods in $8$ model-dataset setups. For
comparison, we also present the performance of the original dense model. As we
can see, the proposed BiP approach finds the best winning tickets (in terms of
the highest accuracy) compared to the baselines across all the pruning setups.
Note that a winning ticket refers to a sparse model with test accuracy that is
at least as good as the dense model (see the marked green dash line). Among
the baseline methods, IMP is the pruning method most competitive to ours.
However, the improvement brought by BiP is significant with respect to the
variance of IMP, except the $60\%$-$80\%$ sparsity regime in (CIFAR-10,
ResNet-20). In the case of (CIFAR-100, ResNet-20) where IMP can not find any
winning tickets (as confirmed by [ma2021sanity]), BiP still manages to find
winning tickets with around $0.6\%$ improvement over the dense model. YY: We
only show BiP and IMP performance on ImageNet due to the computational
resource constraints. In Tab. 1, we summarize the sparest winning tickets
along the pruning trajectory identified by different pruning methods. As we
can see, BiP can identify the winning tickets with higher sparsity levels than
the other methods, except the case of (CIFAR-10, ResNet-20). We also refer
readers to Appendix LABEL:app:_add_results for more experimental results to
justify the high quality of winning tickets identified by BiP.

Fig. 4 demonstrates the structured pruning trajectory on the CIFAR-10/100
datasets. Here we focus on filter pruning, where filter is regarded as a
masking unit in (2). We refer readers to Appendix LABEL:app:_add_results for
channel-wise pruning results. YY: The other experimental results, without
further mentioning, belong to the unstructured pruning regime. In comparison
with Fig. 3, Fig. 4 shows that it becomes much more difficult to find winning
tickets of high-accuracy and high-sparsity in the structured pruning regime,
and the gap among different pruning methods decreases. This is not surprising
since filter pruning imposes much stricter pruning structure constraints than
irregular weight pruning. However, BiP still outperforms all the other
baselines. Most important, it can identify clear winning tickets in the low
sparse regime even when IMP fails.

Figure 5: Time consumption comparison on CIFAR-10 with ResNet-18 to different
pruning ratio $p$.

####  BiP is computationally efficient.

Table 2: Computation complexities of different pruning methods on (CIFAR-10,
ResNet-18). The training epoch numbers of pruning/retraining baselines are
consistent with their official settings or the latest benchmark
implementations. All the evaluations are based on a single Tesla-V100 GPU. SL:
Consider to put in appendix.

Runtime v.s. targeted sparsity |  Training epoch #
---|---
Method $\backslash$ Sparsity | $20\%$ | $59\%$ | $83.2\%$ | $95.6\%$
IMP | 69 min | 276 min | 621 min | 966 min | 160 epoch retrain
Grasp | 89 min | 200 epoch retrain
OMP | 69 min | 160 epoch retrain
Hydra | 115 min |  | 100 epoch prune
---
100 epoch retrain
BiP | 86 min | 100 epoch

In our experiments, another key observation is that BiP yields sparsity-
agnostic runtime complexity while IMP leads to runtime exponential to the
targeted sparsity. Fig. 5 shows the computation cost of different pruning
methods versus different pruning ratios on (CIFAR-10, ResNet-18). For example,
BiP takes 86 minutes (mins) to find the sparsest winning ticket (with $89\%$
sparsity ratio in Tab. 1). This yields $7\times$ less runtime than IMP, which
consumes $620$ mins to find a comparable winning ticket with $87\%$ sparsity.
Compared to the optimization-based pruning baseline Hydra, BiP is also more
efficient as it does not rely on the extra post-optimization model retraining;
see Tab. 2 for a detailed summary of runtime and training epoch numbers
requested by different pruning methods. Further, BiP takes about 1.25 $\times$
more computation time than Grasp and OMP, however, the latter methods lead to
poorer pruning accuracy, e.g., justified by the failure cases of finding
winning tickets in Tab. 1, Fig. 3, and Fig. 4.

#### Rewinding effect on BiP.

Another advantage of BiP is that it is not sensitive to model rewinding.
Recall that rewinding is a strategy used in LTH [frankle2019stabilizing] to
determine what model initialization should be used for re-training a pruned
sparse model. As shown in [ma2021sanity], an IMP-identified winning ticket
could be sensitive to the choice of a rewinding point. Fig. 6 shows the test
accuracy of the BiP-pruned model when it is retrained at different rewinding
epochs, where ‘N/A’ in the x-axis represents the case of no retraining (and
thus no rewinding). As we can see, a carefully-tuned rewinding scheme does not
lead to a significant improvement over BiP without retraining. This suggests
that the subnetworks found by BiP are already of high quality, and justifies
the improved computation efficiency of BiP over Hydra (that requires post-
optimization model retraining), as shown in Fig. 5.

Figure 6: The sensitivity of BiP to rewinding epochs on (CIFAR-10, ResNet-20).
Each curve represents a certain sparsity level and "N/A" in the x-asis
indicates BiP without retraining.

#### Sanity check of BiP on specialized hyperparameters.

Figure 7: The sensitivity of BiP to the lower-level step number $N$ (Left) and
the lower-level regularization parameter $\gamma$ (Right) at different pruning
ratios on (CIFAR-10, ResNet-20).

In Fig.7, we show the sensitivity of BiP to its specialized hyperparameters at
lower-level optimization, including the SGD number ($N$) when customizing
($\boldsymbol{\theta}$-step), and the regularization parameter $\gamma$ in
(2). We refer readers to Fig. LABEL:fig:_ablation_study_full for a more
extensive study. Fig 7(a) shows the test accuracy of BiP-pruned models versus
the choice of $N$ at different pruning ratios. As we can see, the calling of a
larger number of SGD steps at lower-level optimization does not necessarily
improve the overall performance of BiP. This is because in BiP, the
$\boldsymbol{\theta}$-step is initialized by a pre-trained dense model which
does not ask for aggressive weight updating. The best performance of BiP is
typically achieved at $N\leq 3$. We choose $N=1$ throughout the experiments as
it is computationally lightest. Fig 7(b) shows the performance of BiP against
the varying $\gamma$, where ‘$\infty$’ in the x-axis refers to the case of IG-
excluded (3) (by considering $\gamma\to\infty$ in $\mathbf{m}$-step). As we
can see, the choice of $\gamma\in\\{0.5,1\\}$ yields the best pruning accuracy
across all pruning ratios. If $\gamma$ is too small, then the lack of strong
convexity of lower-level optimization would hamper the convergence. If
$\gamma$ is too large, then lower-level optimization would largely departs
from the true model training objective and thus causes a performance drop. In
addition, the accuracy loss encountered at $\gamma=\infty$ also indicates the
necessity of the IG enhancement in BiP.

##  6 Conclusion

We proposed bi-level optimization-enabled pruning method (BiP) to find sparse
networks more effectively and efficiently. To the best of our knowledge, we
are the first to review the model pruning problem from a bi-level point of
view. We theoretically showed that BiP can be solved as easily as the first-
order optimization thanks to the bi-linear nature of the model parameters and
the pruning masks. We conducted extensive experiments to show that BiP can
achieve not only even higher test accuracy than the predominant pruning method
IMP, but also comparable efficiency to the one-shot pruning schemes.

## Checklist

The checklist follows the references. Please read the checklist guidelines
carefully for information on how to answer these questions. For each question,
change the default [TODO] to [Yes] , [No] , or [N/A] . You are strongly
encouraged to include a justification to your answer, either by referencing
the appropriate section of your paper or providing a brief inline description.
For example:

  * •

Did you include the license to the code and datasets? [Yes] See Section
LABEL:gen_inst.

  * •

Did you include the license to the code and datasets? [No] The code and the
data are proprietary.

  * •

Did you include the license to the code and datasets? [N/A]

Please do not modify the questions and only use the provided macros for your
answers. Note that the Checklist section does not count towards the page
limit. In your paper, please delete this instructions block and only keep the
Checklist section heading above along with the questions/answers below.

  1. 1.

For all authors…

    1. (a)

Do the main claims made in the abstract and introduction accurately reflect
the paper’s contributions and scope? [TODO]

    2. (b)

Did you describe the limitations of your work? [TODO]

    3. (c)

Did you discuss any potential negative societal impacts of your work? [TODO]

    4. (d)

Have you read the ethics review guidelines and ensured that your paper
conforms to them? [TODO]

  2. 2.

If you are including theoretical results…

    1. (a)

Did you state the full set of assumptions of all theoretical results? [TODO]

    2. (b)

Did you include complete proofs of all theoretical results? [TODO]

  3. 3.

If you ran experiments…

    1. (a)

Did you include the code, data, and instructions needed to reproduce the main
experimental results (either in the supplemental material or as a URL)? [TODO]

    2. (b)

Did you specify all the training details (e.g., data splits, hyperparameters,
how they were chosen)? [TODO]

    3. (c)

Did you report error bars (e.g., with respect to the random seed after running
experiments multiple times)? [TODO]

    4. (d)

Did you include the total amount of compute and the type of resources used
(e.g., type of GPUs, internal cluster, or cloud provider)? [TODO]

  4. 4.

If you are using existing assets (e.g., code, data, models) or
curating/releasing new assets…

    1. (a)

If your work uses existing assets, did you cite the creators? [TODO]

    2. (b)

Did you mention the license of the assets? [TODO]

    3. (c)

Did you include any new assets either in the supplemental material or as a
URL? [TODO]

    4. (d)

Did you discuss whether and how consent was obtained from people whose data
you’re using/curating? [TODO]

    5. (e)

Did you discuss whether the data you are using/curating contains personally
identifiable information or offensive content? [TODO]

  5. 5.

If you used crowdsourcing or conducted research with human subjects…

    1. (a)

Did you include the full text of instructions given to participants and
screenshots, if applicable? [TODO]

    2. (b)

Did you describe any potential participant risks, with links to Institutional
Review Board (IRB) approvals, if applicable? [TODO]

    3. (c)

Did you include the estimated hourly wage paid to participants and the total
amount spent on participant compensation? [TODO]


# Isoptics of log-aesthetic curves

Ferenc Nagy  University of Debrecen
Doctoral School of Informatics
26 Kassai St.
4028 Debrecen
Hungary

(January 2021)

###### Abstract

The log-aesthetic curve has a significant factor in the field of aesthetic
design to meet the high industrial requirements. It has much exceptional
property and a large number of research papers are published, since its
introduction. It can be parametrized by tangential angle, and it is
autoevolute, meaning that its evolute is also a log-aesthetic curve. In this
paper, we will examine the isoptics of the log-aesthetic curve and determine
if it is autoisoptic as well. We will show that the log-aesthetic curve only
coincides its isoptic in some specific cases. We will also provide a new
formula to calculate the isoptic curve of any tangential angle parametrized
curves.

##  1 Introduction

To design highly aesthetic objects, the use of high-quality curves is very
important. The main property of such curves is the curvature, which should be
monotonically varying. Among them, the log-aesthetic stated as the most
promising for aesthetic design [10] because its curvature and torsion graphs
are both linear [19].

The log-aesthetic curve is originated from Harada et al. [7, 6]. They revealed
that natural aesthetic curves have such a property that their logarithmic
distribution diagram of curvature (LDDC) can be approximated by straight
lines. Based on their work Miura et al. [12, 11] have defined the logarithmic
curvature graph (LCG), an analytical version of the LDDC as the following:
when a curve (given with arc length $s$ and radius of curvature $\rho$) is
subdivided into infinitesimal segments such that $\Delta\rho/\rho$ is
constant, the LCG represents the relationship between $\rho$ and $\Delta s$ in
a double logarithmic graph. The authors formulated the curve whose LCG is
strictly expressed by a straight line (with slope $\alpha$). Such curves
satisfy the following equation:

$\log\rho\frac{\operatorname{d\\!}s}{\operatorname{d\\!}\rho}=\alpha\log\rho+c,$
(1)

where $c$ is a constant. (1) is the fundamental equation of log-aesthetic
curves.

The log-aesthetic curve can also be considered as the generalization of the
clothoid, logarithmic spiral, Nielsen’s spiral, and the circle involute. The
LCG of these curves is a straight line with different slopes that can be drawn
as log-aesthetic curves with different shape parameter $\alpha$ (see Figure
1). The general equations of the log-aesthetic curve will be presented in
Section 2.1.

Figure 1:  The logarithmic curvature graph of several log-aesthetic curves
with different $\alpha$ slopes: clothoid ($\alpha=-1$), Nielsen’s spiral
($\alpha=0$), logarithmic spiral ($\alpha=1$), circle involute ($\alpha=2$),
circle ($\alpha=\raise 1.0625pt\hbox{\set@color$\scriptstyle\pm$}\infty$). The
curves are generated using parametric equations (on the same $[-3\pi,3\pi]$
interval), which parameter is limited (red dots), except in the cases of
$\alpha=1$ and $\alpha=\raise
1.0625pt\hbox{\set@color$\scriptstyle\pm$}\infty$.

Besides the linear curvature, the log-aesthetic curve has another interesting
property. It is autoevolute, that means the evolute of a log-aesthetic curve
with shape parameter $\alpha$ is also a log-aesthetic curve with shape
parameter $-1/(\alpha-2)$ [21].

This characteristic led us to study the isoptic curve of the log-aesthetic
curve. The isoptic curve of a curve is constructed by the involving lines with
a given angle intersect each other at a certain point of the isoptic curve.
Some results of calculating the isoptics of various plane curves will be
presented in Section 2.2. The similar attribute to autoevolute is called the
autoisoptic, meaning that the given curve coincide its isoptic curve. In [13],
there are some autoisoptic curves, among the logarithmic spiral, which is the
log-aesthetic curve with $\alpha=1$. It is also evident that the isoptic of a
circle is also a circle (log-aesthetic curve with $\alpha=\pm\infty$).

In Section 3, we will present a new formula for calculating the isoptic curve
of the log-aesthetic curve (and for all the tangential angle parametrized
curves). We will also prove that for arbitrary $\alpha\in\mathbb{R}$ the log-
aesthetic curve is not autoisoptic.

##  2 Previous results

In this section, we are going to summarize the work of Yoshida and Saito [20].
They analyzed the properties of the log-aesthetic curve and derived a general
formula to draw it. Their related results will be presented in the first
subsection. Besides, in the second part, we will introduce the isoptic curves
and show some results in the topic.

###  2.1 Generating log-aesthetic curves

The general formula is obtained from the relationship between the arc length
and the radius of curvature of the log-aesthetic curve, using a reference
point (at the origin) and placing the following constraints on it: the tangent
vector is $\big{[}1$ $0\big{]}^{T}$, and the radius of curvature is $1$. By
transforming the curve such that the above constraints are satisfied, a point
$P(\theta)$ of the log-aesthetic curve whose tangential angle is $\theta$ can
be expressed on the complex plane as ([20]):

$P(\theta)=\left\\{\begin{array}[]{@{}l@{\quad}l@{}}\int_{0}^{\theta}\mathrm{e}^{(1+{i\mkern
1.0mu})\Lambda\psi}\ d\psi&\text{if }\alpha=1\\\
\int_{0}^{\theta}\big{(}(\alpha-1)\Lambda\psi+1\big{)}^{\frac{1}{\alpha-1}}\mathrm{e}^{{i\mkern
1.0mu}\psi}\ d\psi&\text{otherwise},\end{array}\right.$ (2)

where $\alpha\in\mathbb{R}$ and $\Lambda\in\mathbb{R}^{+}$ are parameters.
$\alpha$ is the slope of the LCG. When $\alpha\neq 1$, all the aesthetic
curves are congruent under similarity transformations without depending on the
value of $\Lambda(\neq 0$). Since similarity transformations preserve angles,
in our examples we focus on the cases of $\Lambda=1$.

The radius of curvature as the function of $\theta$ can be expressed as
([20]):

$\rho(\theta)=\left\\{\begin{array}[]{@{}l@{\quad}l@{}}\mathrm{e}^{\Lambda\theta}&\text{if
}\alpha=1\\\
\big{(}(\alpha-1)\Lambda\theta+1\big{)}^{\frac{1}{\alpha-1}}&\text{otherwise}.\end{array}\right.$
(3)

When $\theta=0$, $\rho=1$. The function increases monotonically when
$\Lambda\neq 0$. In case of $\Lambda=0$, $\rho$ is constant $1$ and the log-
aesthetic curve is a circle. Similarly, when $\alpha=\raise
1.18056pt\hbox{\set@color$\scriptstyle\pm$}\infty$ because taking the limit of
(3) when $\alpha$ approaches $\raise
1.18056pt\hbox{\set@color$\scriptstyle\pm$}\infty$ we get $\rho=1$.

The tangential angle $\theta$ and the arc length $s$ are related by ([20]):

$\theta(s)=\left\\{\begin{array}[]{@{}l@{\quad}l@{}}\frac{1-\mathrm{e}^{-\Lambda
s}}{\Lambda}&\text{if }\alpha=0\\\ \frac{\log(\Lambda s+1)}{\Lambda}&\text{if
}\alpha=1\\\ \frac{(\Lambda\alpha
s+1)^{(1-\frac{1}{\alpha})}-1}{\Lambda(\alpha-1)}&\text{otherwise}.\end{array}\right.$
(4)

Therefore, a point $C(s)$ on the aesthetic curve whose arc length is $s$, can
be defined on the complex plane as ([20]):

$C(s)=\left\\{\begin{array}[]{@{}l@{\quad}l@{}}\int_{0}^{s}\exp({{i\mkern
1.0mu}\frac{1-\mathrm{e}^{-\Lambda u}}{\Lambda}})\ du&\text{if }\alpha=0\\\
\int_{0}^{s}\exp({{i\mkern 1.0mu}\frac{\log(\Lambda u+1)}{\Lambda}})\
du&\text{if }\alpha=1\\\ \int_{0}^{s}\exp({{i\mkern 1.0mu}\frac{(\Lambda\alpha
u+1)^{(1-\frac{1}{\alpha})}-1}{\Lambda(\alpha-1)}})\
du&\text{otherwise}.\end{array}\right.$ (5)

(5) and (2) represent the same curve. Using (4), the radius of curvature can
also be expressed from the arc length $s$ ([20]):

$\rho(s)=\left\\{\begin{array}[]{@{}l@{\quad}l@{}}\mathrm{e}^{\Lambda
s}&\text{if }\alpha=0\\\ \big{(}\Lambda\alpha
s+1\big{)}^{\frac{1}{\alpha}}&\text{otherwise.}\end{array}\right.$ (6)

Since $\rho$ can change from $-\infty$ to $+\infty$, the tangential angle
$\theta$ and arc length $s$ may have upper or lower bound depending on the
value of $\alpha$ because of the negative bases of the fractional exponents
(see Table 1).

Besides the general formula, the authors of [20] also analyzed the properties
of the curve. They called the case of $\Lambda=1$ as the standard form I of
the curve. Let us highlight some important characteristics in this case
regarding the isoptic of the curve, based on the value of $\alpha$:

$\alpha>1$

As $\theta$ approaches $\infty$, the curve spirally diverges towards the point
at $\rho=\infty$ (the arc length to that point is infinite). However, $\theta$
has a lower bound of $\frac{1}{\Lambda(1-\alpha)}$ (see Figure 2(a)–2(d)).

$\alpha<1$

As $\theta$ approaches $-\infty$, the curve spirally converges to the point at
$\rho=0$ (the arc length to that point is infinite). However, $\theta$ has an
upper bound of $\frac{1}{\Lambda(1-\alpha)}$ (see Figure 2(f)–2(l)).

$\alpha=1$

As $\theta$ approaches $-\infty$, the curve spirally converges to the point at
$\rho=0$ (the arc length to that point is finite), and as $\theta$ approaches
$\infty$, the curve spirally diverges towards the point at $\rho=\infty$ (the
arc length to that point is infinite, see Figure 2(e)).

$\alpha=\raise 1.18056pt\hbox{\set@color$\scriptstyle\pm$}\infty$

The log-aesthetic curve is a circle centered at $\big{[}0$ $1\big{]}^{T}$.

| Tangential angle ($\theta$) |  | Arc length ($s$)
---|---|---|---
| $\alpha<1$ | $\alpha=1$ | $\alpha>1$ |  | $\alpha<0$ | $\alpha=0$ | $\alpha>0$
Upper bound: | $\frac{1}{\Lambda(1-\alpha)}$ | - | - |  | $-\frac{1}{\Lambda\alpha}$ | - | -
Lower bound: | - | - | $\frac{1}{\Lambda(1-\alpha)}$ |  | - | - | $-\frac{1}{\Lambda\alpha}$
Table 1: Upper and lower bound of $\theta$ and $s$.

(a) $\alpha=3$

(b) $\alpha=2.5$

(c) $\alpha=2$,
Circle involute

(d) $\alpha=1.5$

(e) $\alpha=1$,
Logarithmic spiral

(f) $\alpha=0.5$

(g) $\alpha=0$,
Nielsen’s spiral

(h) $\alpha=-0.5$

(i) $\alpha=-1$,
Clothoid

(j) $\alpha=-1.5$

(k) $\alpha=-2$

(l) $\alpha=-2.5$

Figure 2: Log-aesthetic curves with different shape parameter $\alpha$ and
$\Lambda=1$. The red points are the current bounds of $\theta$.

###  2.2 Isoptic curves

###### Definition 1

For a given plane curve the isoptic curve is the locus of points from where
the given curve can be seen under a predefined angle ${(0<\gamma<\pi)}$[18].
The isoptic curve with right angle called orthoptic curve.

For some planar, second-degree curves, the isoptic can be constructed simply
based on the definition (see Figure 3, or for an overview in [18]). The
authors of [14] calculated isoptic curves of conic sections. There are also
interesting results of curves whose isoptic is a circle in [16, 9], or whose
isoptic is an ellipse in [17]. The work of Cieślak et al. is also remarkable
[2, 3], in which, the authors define the isoptics of closed convex curves
using support function. Besides, there are also important results of isoptics
in non-Euclidean geometry, e.g. in [5] or in [4].

Figure 3: Isoptic curves (green) of a line segment and a parabola with
${\gamma=\frac{3\pi}{4}}\textrm{, }{\frac{2\pi}{3}}\textrm{,
}{\frac{\pi}{2}}\textrm{, and }{\frac{\pi}{3}}$. The orthoptic of the line
segment is a circle, and the ortoptic of the parabola is its directrix.

It should be noted that for higher degree curves the direct computation of the
isoptic based on the definition can be extremely difficult. Therefore, Kunkli
et al. defined a new calculation for free-form curves [8]. The authors
determine the isoptic curve of convex Bézier-curves as the envelope of
envelopes of families of isoptic circles over the chords of the curve.
Although, the calculations can easily be applied for other curve types, such
as B-splines or NURBS, for complex or higher degree curves the isoptic usually
can only be approximated.

##  3 Isoptics of log-aesthetic curves

Based on Definition 1 and using the tangential angle parametrized (2) of the
log-aesthetic curve the isoptic curve can be constructed as the intersection
of the appropriate involving lines. For an arbitrary point $P(\theta)$ of the
log-aesthetic curve the corresponding curve point, which tangent line
intersects the involving line of $P(\theta)$ under a given angle $\gamma$ is
calculated as $P\left(\theta+(\pi-\gamma)\right)$. The tangent with unit
length of the point $\theta$ is: $[\cos\theta,\sin\theta]^{T}$. It is also
true for all the tangential angle parametrized planar $r(t)$ curves [1]:

$\displaystyle\frac{r^{\prime}(t)}{\left|r^{\prime}(t)\right|}$
$\displaystyle=\begin{bmatrix}\cos t\\\ \sin t\end{bmatrix}.$ (7)

Let us define the isoptic with angle $\delta=\pi-\gamma$. Thus, the parametric
equations of the appropriate tangent lines are

$\displaystyle\begin{split}P(t_{\theta})&=P(\theta)+t_{\theta}\begin{bmatrix}\cos\theta\\\
\sin\theta\end{bmatrix}\\\
P(t_{\theta+\delta})&=P(\theta+\delta)+t_{\theta+\delta}\begin{bmatrix}\cos(\theta+\delta)\\\
\sin(\theta+\delta)\end{bmatrix}.\end{split}$ (8)

Solving it for $t_{\theta}$, we obtain:

$t_{\theta}=\csc(\delta)\left(\vec{V}^{\delta}_{x}(\theta)\sin(\theta+\delta)-\vec{V}^{\delta}_{y}(\theta)\cos(\theta+\delta)\right),$
(9)

where $\vec{V}^{\delta}_{x}(\theta)$ and $\vec{V}^{\delta}_{y}(\theta)$ are
the coordinates of the vector
$\vec{V}^{\delta}(\theta)=\overrightarrow{P(\theta)P(\theta+\delta)}$, which
explicit equation can be written (from (2)) on the complex plane as

$\vec{V}^{\delta}(\theta)=\left\\{\begin{array}[]{@{}l@{\quad}l@{}}\int_{\theta}^{\theta+\delta}\mathrm{e}^{(1+{i\mkern
1.0mu})\Lambda\psi}\ d\psi&\text{if }\alpha=1\\\
\int_{\theta}^{\theta+\delta}\big{(}(\alpha-1)\Lambda\psi+1\big{)}^{\frac{1}{\alpha-1}}\mathrm{e}^{{i\mkern
1.0mu}\psi}\ d\psi&\text{otherwise}.\end{array}\right.$ (10)

Therefore, the isoptic curve of the log-aesthetic (and all the tangential
angle parametrized) curve(s) with a given $\delta$ ($=\pi-\gamma$) angle can
be calculated as (see Figure 4):

$I_{\delta}(\theta)=P(\theta)+\csc(\delta)\left(\vec{V}^{\delta}_{x}(\theta)\sin(\theta+\delta)-\vec{V}^{\delta}_{y}(\theta)\cos(\theta+\delta)\right)\begin{bmatrix}\cos\theta\\\
\sin\theta\end{bmatrix}$ (11)

The following modifications giving us more geometric meaning. The (9) can be
altered using the harmonic addition theorem [15] as

$t_{\theta}=\csc(\delta)\operatorname{sign}\left(\vec{V}^{\delta}_{x}(\theta)\right)\mathinner{\\!\left\lVert\vec{V}^{\delta}(\theta)\right\rVert}\sin{\left(\theta+\delta+\arctan{\left(\frac{-\vec{V}^{\delta}_{y}(\theta)}{\vec{V}^{\delta}_{x}(\theta)}\right)}\right)}.$
(12)

Changing the $\arctan$ function to the well-known $\operatorname{arctan2}$,
that calculates the tangential angle of vector $\vec{V}^{\delta}(\theta)$
(a.k.a. the angle between the positive x axis and the ray to the point
$\big{[}\vec{V}^{\delta}_{x}(\theta),\vec{V}^{\delta}_{y}(\theta)\big{]}^{T}$
and returns the value in the interval of $[0,\pi]$) the
$\operatorname{sign}\left(\vec{V}^{\delta}_{x}(\theta)\right)$ can be omitted:

$t_{\theta}=\csc(\delta)\mathinner{\\!\left\lVert\vec{V}^{\delta}(\theta)\right\rVert}\sin{\left(\theta+\delta-\operatorname{arctan2}{\left(\vec{V}^{\delta}_{y}(\theta),\vec{V}^{\delta}_{x}(\theta)\right)}\right)}.$
(13)

It means that the isoptic curve of the tangential angle parametrized curves
can also be given as

$I_{\delta}(\theta)=P(\theta)+\csc(\delta)\mathinner{\\!\left\lVert\vec{V}^{\delta}(\theta)\right\rVert}\sin{\left(\theta+\delta-\operatorname{arctan2}{\left(\vec{V}^{\delta}_{y}(\theta),\vec{V}^{\delta}_{x}(\theta)\right)}\right)}\begin{bmatrix}\cos\theta\\\
\sin\theta\end{bmatrix}$ (14)

(a) $\alpha=2$
Circle involute

(b) $\alpha=1$
Logarithmic spiral

(c) $\alpha=0$
Nielsen’s spiral

(d) $\alpha=-1$
Clothoid

Figure 4: Log-aesthetic curves with $\Lambda=1$ (blue) and their isoptics with
given $\delta=\frac{\pi}{3}$ (green), with the point of bounds (red and
orange).

To analyze the autoisoptic property, let us see the logarithmic curvature
graph of the log-aesthetic curve, parametrized by $\theta$:

$P_{LCG}(\theta)=\begin{bmatrix}\log\rho(\theta)\\\
\log\frac{\rho(\theta)}{\mid\rho(s)\mid}\end{bmatrix},$ (15)

By substituting the appropriate $\rho$ functions, after some simplifications,
it is:

$P_{LCG}(\theta)=\left\\{\begin{array}[]{@{}l@{\quad}l@{}}\begin{bmatrix}\Lambda\\\
\theta\Lambda-\log(\Lambda)\end{bmatrix}&\alpha=1\\\\[10.00002pt]
\begin{bmatrix}\frac{\log((\alpha-1)\theta\Lambda+1)}{\alpha-1}\\\
\frac{\alpha\log((\alpha-1)\theta\Lambda+1)}{\alpha-1}-\log(\Lambda)\end{bmatrix}&\text{otherwise.}\end{array}\right.$
(16)

The (16) on the LCG (which axes are $\log\rho$ and $\log\frac{\Delta
s}{\Delta\rho/\rho}$) defines a straight line with slope $\alpha$ (see Figure
1). The slope (which is constant for log-aesthetic curves) can be calculated
as the following (by taking into account the Table 1 limits of $\theta$):

$\hat{\alpha}(\theta)=\left\\{\begin{array}[]{@{}l@{\quad}l@{}}\frac{P_{LCG}(\theta-\phi)_{y}-P_{LCG}(\theta)_{y}}{P_{LCG}(\theta-\phi)_{x}-P_{LCG}(\theta)_{x}}&\alpha<1\\\\[10.00002pt]
\frac{P_{LCG}(\theta+\phi)_{y}-P_{LCG}(\theta)_{y}}{P_{LCG}(\theta+\phi)_{x}-P_{LCG}(\theta)_{x}}&\alpha\geq
1,\end{array}\right.$ (17)

where $\phi\in\\{\mathbb{R}\setminus{0}\\}$ and $P_{LCG}(\theta)_{x}$ and
$P_{LCG}(\theta)_{y}$ is the coordinates of the point on the LCG defined by
$\theta$.

Unfortunately, the exact logaritmic curvature graph of the isoptic curve,
using the $\theta$ parametrized (11) or (14) cannot be determined for
arbitrary $\alpha\in\ \mathbb{R}$ (not even using computer algebra softwares).
Solely, in the cases when the log-aesthetic curve is the logarithmic spiral
($\alpha=1$) or when the (2) can be represented in terms of trigonometric
functions ($\alpha=2,\frac{3}{2},\frac{4}{3},\dots,\frac{k+1}{k}$,
$k\in\mathbb{N}^{*}$) or using Fresnel integrals ($\alpha=-1,3$) [22].

Let us see the case of $\alpha=1$. The LCG equation of the log-aestheic curve
is:

$P^{\alpha=1}_{LCG}(\theta)=\begin{bmatrix}\theta\Lambda\\\
\theta\Lambda-\log(\Lambda)\end{bmatrix},$ (18)

and the LCG formula of its isoptic is

$I^{\alpha=1}_{LCG}(\theta)=\begin{bmatrix}\frac{1}{2}\log\left(\frac{e^{2\delta\Lambda}-2e^{\delta\Lambda}\cos(\delta)+1}{\Lambda^{2}+1}\right)+\log(\csc(\delta))+\theta\Lambda\\\
\log\left(\frac{\csc(\delta)\sqrt{\frac{e^{2\delta\Lambda}-2e^{\delta\Lambda}\cos(\delta)+1}{\Lambda^{2}+1}}}{\Lambda}\right)+\theta\Lambda\end{bmatrix}.$
(19)

By calculating the slope $\hat{\alpha}(\theta)$ of the isoptic, using a
similar formula as (17), we see that it is constant $1$ for all
$\theta\in\mathbb{R}$, meaning that, it is also a log-aesthetic curve with
shape parameter 1. It proves its autoisoptic property (see Figure 5).

Figure 5: Logarithmic curvature graph of the log-aesthetic curve with
$\alpha=1$, $\Lambda=1$ (blue) and its isoptic with $\delta=\frac{2}{3}\pi$
(green, dashed).

For another example, let us see the LCG of the log-aesthetic curve with
$\alpha=2$ (circle involute), that is

$P^{\alpha=2}_{LCG}(\theta)=\begin{bmatrix}\log(\theta\Lambda+1)\\\
2\log(\theta\Lambda+1)-\log(\Lambda)\end{bmatrix},$ (20)

and the LCG of the isoptic with $\delta=\frac{2\pi}{3}$ and $\Lambda=1$ is
(note: it is not the general formula with arbitrary $\Lambda>0$ and
${0<\delta<\pi}$ because it is too long to appear in print):

$I^{\alpha=2,\Lambda=1}_{LCG,\delta=2\pi/3}(\theta)=\begin{bmatrix}\log\left(\frac{2\left(-6\pi\left(-3\theta+\sqrt{3}-3\right)+27\theta(\theta+2)+4\pi^{2}+54\right)^{3/2}}{27\sqrt{3}\left|3\theta(\theta+2)+\pi\left(2\theta-\sqrt{3}+2\right)+\frac{4\pi^{2}}{9}+9\right|}\right)\\\\[25.00003pt]
\log\left(\frac{2\left(-6\pi\left(-3\theta+\sqrt{3}-3\right)+27\theta(\theta+2)+4\pi^{2}+54\right)^{3/2}}{27\sqrt{3}}\right.+\\\
+\left.\frac{\left|\frac{-9\pi\left(-2\theta+\sqrt{3}-2\right)+27\theta(\theta+2)+4\pi^{2}+81}{(3\theta+\pi+3)\left(3\pi\left(6\theta-5\sqrt{3}+6\right)+27(\theta(\theta+2)+5)+4\pi^{2}\right)}\right|}{27\sqrt{3}}\right)\end{bmatrix}.$
(21)

By calculating its $\hat{\alpha}(\theta)$ slope, similarly as (17) with
different $\theta$ values (and with $\phi=\pi$), we see:

$\begin{split}\hat{\alpha}(-1)&\approx 1.86204\text{ (lower bound of
}\theta\text{)}\\\ \hat{\alpha}(0)&\approx 1.96486\\\
\hat{\alpha}(\pi)&\approx 1.99523\\\ \hat{\alpha}(3\pi)&\approx 1.99899\\\
\hat{\alpha}(10\pi)&\approx 1.99987\\\ \hat{\alpha}(100\pi)&\approx
1.99999845\\\ \lim_{\theta\to\infty}\hat{\alpha}(\theta)&=2.\end{split}$ (22)

Therefore, it can be claimed that the log-aesthetic curve with $\alpha=2$ does
not coincide with its isoptic, meaning that the circle involute is not
autoisoptic (see Figure 6). However, the limit of $\hat{\alpha}(\theta)$ is 2,
when $\theta$ approaches $\infty$.

Figure 6: Logarithmic curvature graph of the log-aesthetic curve with
$\alpha=2$, $\Lambda=1$ (blue and the dashed, that is a transformed) and its
isoptic with $\delta=\frac{2}{3}\pi$ (green).

We receive similar results by calculating and drawing the different
logarithmic curvature graphs of the isoptic curve in the other mentioned
$\alpha>1$ cases. The slope of the isoptic is not constant, but the limit of
the slope as the tangential angle approaches infinity is always equal to the
specified $\alpha$ parameter. We can draw the following conclusion: the log-
aesthetic curve is not autoisoptic in all the cases of $\alpha>1$.

Let us see another example, the log-aesthetic curve with $\alpha=-1$
(clothoid). In this case, the general equation of the log-aesthetic curve is
defined using Fresnel integrals [22]. The LCG of it can be obtained as

$P^{\alpha=-1}_{LCG}(\theta)=\begin{bmatrix}-\frac{1}{2}\log(1-2\theta\Lambda)\\\
\frac{1}{2}\log(1-2\theta\Lambda)-\log(\Lambda)\end{bmatrix}.$ (23)

However, the LCG of its isoptic can be obtained with a fixed $\Lambda$ and
$\delta$ values using computer algebra systems, it is too long to appear in
print. We calculated it for $\Lambda=1$ and $\delta=\frac{2}{3}\pi$ and
computed the $\hat{\alpha}(\theta)$ slope of it with different $\theta$ values
and with $\phi=\pi$, similarly as in (17). From the following cases

$\begin{split}\hat{\alpha}(\frac{1}{2}-\delta)&\approx-2.23102\text{ (upper
bound of }\theta\text{ for the isoptic)}\\\
\hat{\alpha}(-\pi)&\approx-1.27031\\\ \hat{\alpha}(-2\pi)&\approx-1.0705\\\
\hat{\alpha}(-5\pi)&\approx-1.01251\\\ \hat{\alpha}(-10\pi)&\approx-1.00329\\\
\hat{\alpha}(-100\pi)&\approx-1.00003\\\
\lim_{\theta\to-\infty}\hat{\alpha}(\theta)&=-1,\end{split}$ (24)

we see that the log-aesthetic curve with $\alpha=-1$ does not coincide with
its isoptic, meaning that the clothoid is not autoisoptic (see Figure 7).
However, the limit of $\hat{\alpha}(\theta)$ is -1, when $\theta$ approaches
$-\infty$. Note, the slope $\hat{\alpha}(\theta)$ of the LGC of the isoptic in
this case differs more from $\alpha$ around the upper bound rather than in the
case of $\alpha=2$. The reason can be that fact the log aesthetic curve has an
inflection point at the upper bound of $\theta$ in case of $\alpha<1$.

Figure 7: Logarithmic curvature graph of the log-aesthetic curve with
$\alpha=-1$, $\Lambda=1$ (blue and the dashed, that is a transformed) and its
isoptic with $\delta=\frac{2}{3}\pi$ (green).

##  4 Conclusion

We have defined the isoptics of the log-aesthetic curve in the two-dimensional
plane. We have also given an exact formula to obtain the isoptics of all the
tangential angle parametrized plane curves.

We have examined three instances of the log-aesthetic curve based on its shape
parameter. In the first case, when it is the logarithmic spiral ($\alpha=1$),
we have shown its autoisoptic property. In the second case, when it is the
circle involute ($\alpha=2$) and in the third case, when it is the clothoid
($\alpha=-1$), we have proven that they are not autoisoptic curves. It means
that the log-aesthetic curve with arbitrary $\alpha\in\mathbb{R}$ is not
autoisoptic.

Moreover, by analyzing (13), we believe that the log-aesthetic curve coincides
with its isoptic only in the cases of $\alpha=1$ and $\alpha=\raise
1.18056pt\hbox{\set@color$\scriptstyle\pm$}\infty$ because in this cases the
tangential angle $\theta$ and the tangential angle of the vector
$\vec{V}^{\delta}(\theta)$ are in a relation of proportionality
($\theta\propto\operatorname{arctan2}{\left(\vec{V}^{\delta}_{y}(\theta),\vec{V}^{\delta}_{x}(\theta)\right)}$),
such that the
$\sin{\left(\theta+\delta-\operatorname{arctan2}{\left(\vec{V}^{\delta}_{y}(\theta),\vec{V}^{\delta}_{x}(\theta)\right)}\right)}$
is constant and the
$\mathinner{\\!\left\lVert\vec{V}^{\delta}(\theta)\right\rVert}$ equals to the
radius of curvature of the log-aesthetic curve with a constant multiplier. In
the other cases, when $-\infty<\alpha<1$ or $1<\alpha<\infty$ the ratio of the
tangential angle $\theta$ and the tangential angle of the vector
$\vec{V}^{\delta}(\theta)$ is not constant, therefore the log-aesthetic curve
does not fulfill these properties. We assume that this is the reason why the
log-aesthetic curve is not autoisoptic in these cases.

## References

  * [1] Elsa Abbena, Simon Salamon, and Alfred Gray.  Modern differential geometry of curves and surfaces with Mathematica.  CRC press, 1998.
  * [2] Waldemar Cieślak, Andrzej Miernowski, and Witold Mozgawa.  Isoptics of a closed strictly convex curve.  In Global Differential Geometry and Global Analysis, pages 28–35. Springer Berlin Heidelberg, 1991.
  * [3] Waldemar Cieślak, Andrzej Miernowski, and Witold Mozgawa.  Isoptics of a closed strictly convex curve.-ii.  Rendiconti del Seminario Matematico della Università di Padova, 96:37–49, 1996.
  * [4] Géza Csima and Szirmai Jenő.  Isoptic curves of conic sections in constant curvature geometries.  Mathematical Communications, 19(2):277–290, 2014.
  * [5] Géza Csima and Jenő Szirmai.  On the isoptic hypersurfaces in the n-dimensional euclidean space.  KoG (Scientific and professional journal of Croatian Society for Geometry and Graphics), 17(17.):53–57, 2013.
  * [6] T. Harada, F. Yoshimoto, and M. Moriyama.  An aesthetic curve in the field of industrial design.  In Proceedings 1999 IEEE Symposium on Visual Languages, pages 38–47. IEEE, 1999.
  * [7] Toshinobu Harada.  Study of quantitative analysis of the characteristics of a curve.  Forma, 12(1):55–63, 1997.
  * [8] Roland Kunkli, Ildikó Papp, and Miklós Hoffmann.  Isoptics of bézier curves.  Computer Aided Geometric Design, 30(1):78–84, 2013.
  * [9] Árpád Kurusa.  Is a convex plane body determined by an isoptic?  Beiträge zur Algebra und Geometrie/Contributions to Algebra and Geometry, 53(1):281–294, 2012.
  * [10] Raph Levien and Carlo H. Séquin.  Interpolating splines: which is the fairest of them all?  Computer-Aided Design and Applications, 6(1):91–102, 2009.
  * [11] Kenjiro T. Miura.  A general equation of aesthetic curves and its self-affinity.  Computer-Aided Design and Applications, 3(1-4):457–464, 2006.
  * [12] Kenjiro T Miura, Junji Sone, Atsushi Yamashita, and Toru Kaneko.  Derivation of a general formula of aesthetic curves.  In Proceedings of the Eighth International Conference on Humans and Computers (HC2005), pages 166–171, 2005.
  * [13] Boris Odehnal.  Examples of autoisoptic curves.  In Luigi Cocchiarella, editor, ICGG 2018 - Proceedings of the 18th International Conference on Geometry and Graphics, pages 350–361, Cham, 2019. Springer International Publishing.
  * [14] Paul Siebeck.  Ueber eine gattung von curven vierten grades, welche mit den elliptischen functionen zusammenhängen.  Journal für die reine und angewandte Mathematik, 1860(57):359–370, 1860.
  * [15] Eric W. Weisstein.  Harmonic addition theorem. From MathWorld—A Wolfram Web Resource.  Visited: 12/01/2021.
  * [16] Walter Wunderlich.  Kurven mit isoptischem kreis.  Aequationes mathematicae, 6(1):71–81, 1971.
  * [17] Walter Wunderlich.  Kurven mit isoptischer ellipse.  Monatshefte für Mathematik, 75(4):346–362, 1972.
  * [18] R. C. Yates.  A Handbook on Curves and their Properties.  Ann Arbor, J.W. Edwards, 1947.
  * [19] Norimasa Yoshida, Ryo Fukuda, and Takafumi Saito.  Log-aesthetic space curve segments.  In 2009 SIAM/ACM Joint Conference on Geometric and Physical Modeling, pages 35–46, 2009.
  * [20] Norimasa Yoshida and Takafumi Saito.  Interactive aesthetic curve segments.  The Visual Computer, 22(9-11):896–905, 2006.
  * [21] Norimasa Yoshida and Takafumi Saito.  The evolutes of log-aesthetic planar curves and the drawable boundaries of the curve segments.  Computer-Aided Design and Applications, 9(5):721–731, 2012.
  * [22] Rushan Ziatdinov, Norimasa Yoshida, and Tae wan Kim.  Analytic parametric equations of log-aesthetic curves in terms of incomplete gamma functions.  Computer Aided Geometric Design, 29(2):129–140, 2012.