# SHOT: Suppressing the Hessian along the Optimization Trajectory for
Gradient-Based Meta-Learning

JunHoo Lee, Jayeon Yoo, and Nojun Kwak
Seoul National University
{mrjunoo, jayeon.yoo<EMAIL_ADDRESS>Corresponding author.

###### Abstract

In this paper, we hypothesize that gradient-based meta-learning (GBML)
implicitly suppresses the Hessian along the optimization trajectory in the
inner loop. Based on this hypothesis, we introduce an algorithm called SHOT
(Suppressing the Hessian along the Optimization Trajectory) that minimizes the
distance between the parameters of the target and reference models to suppress
the Hessian in the inner loop. Despite dealing with high-order terms, SHOT
does not increase the computational complexity of the baseline model much. It
is agnostic to both the algorithm and architecture used in GBML, making it
highly versatile and applicable to any GBML baseline. To validate the
effectiveness of SHOT, we conduct empirical tests on standard few-shot
learning tasks and qualitatively analyze its dynamics. We confirm our
hypothesis empirically and demonstrate that SHOT outperforms the corresponding
baseline. Code is available at: https://github.com/JunHoo-Lee/SHOT

##  1 Introduction

With the advent of deep learning revolution in the 2010s, machine learning has
widened its application to many areas and outperformed humans in some areas.
Nevertheless, machines are still far behind humans in the ability to learn by
themselves, in that humans can learn concepts of new data with only a few
samples, while machines cannot. Meta-Learning deals with this problem of
learning to learn. An approach to address this issue is the metric-based
methods Chen and He (2021); Snell et al. (2017); Bromley et al. (1993); Grill
et al. (2020); Chen et al. (2020); Sung et al. (2018); Radford et al. (2021),
which aims to learn ‘good‘ kernels in order to project given data into a well-
defined feature space. Another popular line of research is optimization-based
methods Finn et al. (2017); Ravi and Larochelle (2017); Rusu et al. (2018a),
so-called Gradient-Based Meta-Learning (GBML), which aims to achieve the goal
of meta-learning by gradient descent.

As the most representative GBML methods applicable to any model trained with
the gradient descent process, model-agnostic meta-learning (MAML) Finn et al.
(2017) and its variations Raghu et al. (2019); Oh et al. (2020); Rusu et al.
(2018b) exploit nested loops to find a good meta-initialization point from
which new tasks can be quickly optimized. To do so, it divides the problem
into two loops: the inner loop (task-specific loop) and the outer loop (meta
loop). The former tests the meta-learning property of fast learning, and the
latter moves the meta-initialization point by observing the dynamics of the
inner loop. The general training is performed by alternating the two loops. In
this paper, by focusing on the inner loop of GBML methods, we propose a new
loss term for the outer loop optimization.

To properly evaluate the meta-learning property in the inner loop, GBML
samples a new problem for each inner loop. Fig. 1 shows how the problem is
formulated. For every inner loop, the task is sampled from the task
distribution. The classes constituting each task are randomly sampled, and the
configuration order (class order) is also randomized. This means that the
model tackles an entirely new task at the beginning of each inner loop. Since
the task is unknown and the model has to solve the problem within a few
gradient steps (some extreme algorithms Raghu et al. (2019); Oh et al. (2020)
solve the problem in one gradient step), the model typically exploits large
learning rates to adapt quickly. For example, MAML Finn et al. (2017) and ANIL
Raghu et al. (2019) use significantly larger learning rates of 0.4 and 0.5
respectively, for the inner loop compared to that of the outer loop
($10^{-3}$). These rates are much larger than those used in modern deep
learning techniques such as ViT Dosovitskiy et al. (2020), which typically use
rates in the order of $10^{-6}$. However, using such large learning rates can
lead to instability due to high-order terms like the Hessian LeCun et al.
(2015). The issue of instability caused by high-order terms is particularly
evident during the initial training phase. To mitigate this problem, some
deep-learning methodologies employ a heuristic that utilizes a smaller
learning rate at the beginning Dosovitskiy et al. (2020); He et al. (2018).
Given that GBML goes through only several optimizations in the entire inner
loop, GBML’s inner loop corresponds to the initial training phase. Therefore,
we anticipate that the problem caused by high-order terms may be severe with a
large learning rate in each inner loop. However, it is noteworthy that GBML
appears to perform well despite these potential issues.

In this paper, we investigate the inner loop of GBML and its relationship with
the implicit risk of using a large learning rate with a few optimization steps
in the inner loop. We observe that the gradient of a GBML model remains
relatively constant within an inner loop, indicating that the Hessian along
the optimization trajectory has a minimal impact on the model’s dynamics as
learning progresses. We hypothesize that this phenomenon is a key property
that enables GBML to mitigate the risk of using a large learning rate in the
inner loop.

Our hypothesis that GBML implicitly suppresses the Hessian along the
optimization trajectory suggests a desirable property for a GBML algorithm. By
explicitly enforcing this property instead of relying on implicit suppression,
we anticipate that the model will achieve faster convergence and superior
performance. However, enforcing this property is not trivial as it involves
dealing with high-order terms, which naturally require a significant amount of
computation. In this paper, we propose an algorithm called Suppressing the
Hessian along the Optimization Trajectory (SHOT) that suppresses the Hessian
in the inner loop by minimizing the distance of the target model from a
reference model which is less influenced by the Hessian. Our algorithm can
enforce the Hessian to have desirable properties along the inner loop without
much increasing the computational complexity of a baseline model. More
specifically, SHOT only requires one additional forward pass without any
additional backward pass in general cases.

Our algorithm, SHOT, is agnostic to both the algorithm and architecture used
in GBML, making it highly versatile and applicable to any GBML baseline. To
validate the effectiveness of SHOT, we conducted empirical tests on standard
few-shot learning tasks, including miniImagenet, tieredImagenet, Cars, and CUB
Vinyals et al. (2016); Ren et al. (2018); Krause et al. (2013); Welinder et
al. (2010). Furthermore, we tested SHOT’s cross-domain ability by evaluating
its performance on different benchmarks. Also, SHOT can be applied to one-step
or Hessian-free baselines, acting as a regularizer. Our results demonstrate
that SHOT can significantly boost GBML in an algorithm- and architecture-
independent manner. Finally, we analyzed the dynamics of SHOT and confirmed
that it behaves as expected.

Figure 1:  For the $i$-th inner loop in meta-learning, the task, $\tau_{i}$,
is sampled from the task distribution: $N$ (e.g. 3) classes constituting each
task are randomly sampled with a random configuration order, meaning that the
model does not have any information before it retrieves an unseen task.

##  2 Related Works

GBML, or optimization-based meta-learning, is a type of methodology which
wants to mimic people’s rapid learning skill. Thus, a model should learn with
a small number of samples and gradient steps. To achieve this property, Ravi
and Larochelle (2017) tried to save meta-data to recurrent networks. Recently,
MAML-based methods Finn et al. (2017); Rajeswaran et al. (2019); Rusu et al.
(2018b); Baik et al. (2020); Bernacchia (2021) are typically considered as a
synonym of GBML consisting of nested loops: outer loop (meta-loop) and inner
loop (task-specific loop). They sample few-shot classification tasks to test
the meta-learning property in the inner loop and update the meta-
initialization parameters in the outer loop with SGD. To learn rapidly, they
exploit a more-than-thousand times larger learning rate in the inner loop
compared to modern deep learning settings such as ViT Dosovitskiy et al.
(2020) and CLIP Radford et al. (2021). Although some researches have proved
their convergence property Fallah et al. (2020); Wang et al. (2020), they set
infinitesimal size of learning rate to prove convergence. Fallah et al. (2020)
set the learning rate considering the Hessian i.e., setting the learning rate
small enough so that the higher-order terms do not affect the convergence and
Wang et al. (2020) proved the convergence with the assumption that the
learning rate goes to zero. Although some lines of research have shown
algorithms dealing with learning rates Baik et al. (2020); Bernacchia (2021),
they did not explain the success of a large learning rate. Baik et al. (2020)
used an additional model to predict a proper task-specific learning rate. ANIL
(Almost No Inner Loop) Raghu et al. (2019) and BOIL (Body Only Inner Loop) Oh
et al. (2020) tried to explain GBML with the dynamics of the encoder’s output
features. Although both used the same feature-level perspective, they argued
in exactly different ways. ANIL Raghu et al. (2019) argued feature reuse is an
essential component, while BOIL Oh et al. (2020) said that feature adaption is
crucial. Both algorithms can be considered as a variant of preconditioning
meta-learning since they freeze most layers in the inner loop, thereby frozen
layers can be considered as warp-parameters Flennerhag et al. (2019). Our SHOT
can be viewed as a variant of preconditioning, as we hypothesize that there
exists a suitable condition of meta-parameter that can adapt well to a given
meta-learning task. However, unlike previous work, we define the condition of
the meta-parameter based on the curvature of the loss surface. Additionally,
our algorithm does not require any specific architecture such as warp-
parameters, making it applicable to any architecture and any GBML algorithm.
Finally, there exists some paper which delt with curvature, those all required
some special architecture Simon et al. (2020), or explicit calculation of the
Hessian Hiller et al. (2022); Park and Oliva (2019). In contrast, we propose
much weaker condition and does not require specific architecture. We only
suppress the Hessian along the inner loop. In addition, we discuss in the
Appendix the possibility of integrating the hypotheses proposed in ANIL Raghu
et al. (2019) and BOIL Oh et al. (2020) into our SHOT framework.

Optimization methods based on second-order or first-order derivatives that
take curvature into account is a vast area of research. Although they incur a
higher computational cost than vanilla SGD, they have been demonstrated to
exhibit superior convergence properties. Several studies have investigated the
curvature of the loss surface while maintaining the first-order computational
costs. Natural gradient descent Bonnabel (2013), a type of Riemannian gradient
descent employing the Fisher information metric as a measure, is coordinate-
invariant, and exhibits better convergence properties than vanilla gradient
descent. Given that GBML utilizes relatively fewer samples, incorporating
second-order terms is computationally more feasible. As a result, there is a
growing body of research that employs natural gradient descent in the inner
loop. Additionally, the SAM optimizer Foret et al. (2020) implicitly leverages
curvature by iteratively searching for flat minima during training, leading to
better generalization properties than SGD. This means they find a point where
the Hessian is close to zero around a local minimum. Our method can be
interpreted as the variation of Knowledge Distillation Hinton et al. (2015).
We reduced the distance between models, slow model and fast model. However, we
do not need any pretrained model. Those models share same meta-initialization
parameter.

As our goal is to find a proper meta-initialization point that is not affected
by the Hessian along the inner loop, our proposed algorithm, SHOT, also deals
with the Hessian of the loss surface. However, unlike SAM optimizer Foret et
al. (2020), we do not seek a point where the Hessian is zero in any direction.
We rather seek a point where the Hessian does not take effect only along the
optimization trajectory in the inner loop. By using this restriction, we give
a looser constraint. Also, while the SAM optimizer has a nested loop to find a
proper point, SHOT does not need any iteration in a single optimization step.
Also, although natural gradient descent regards curvature like our algorithm,
our SHOT does not calculate the exact Hessian during inference time as we
precondition the Hessian and only consider the Hessian’s effect along the
optimization trajectory. In this paper, we aim to shed light on how GBML’s
inner loop works. Through an analysis of the curvature of the loss surface in
the inner loop, we hypothesize that the primary role of the outer loop is to
implicitly regulate the curvature of the inner loop’s loss surface. To
reinforce this property by explicitly incorporating a corresponding loss term
in the inner loop, we propose an algorithm called SHOT, which takes into
account the curvature of the loss surface in the inner loop and aims to
enhance its generalization properties and dynamics by controlling the
curvature itself. We provide empirical evidence supporting our hypothesis
through the success of SHOT. Additionally, we offer a new interpretation of
GBML as a variant of MBML in the Appendix.

##  3 Preliminaries: GBML

Figure 2:  This figure provides the concept of the proposed SHOT.
Specifically, a reference model ($\hat{\theta}^{\star}$) is initialized with
the same support set and starting point as our model (${\theta}^{\star}$), but
undergoes significantly more optimization steps. The outer loop then minimizes
the distance between the two models.

We consider the few-shot classification task of $N$-way $K$-shot, which
requires learning $K$ samples per class to distinguish $N$ classes in each
inner loop. In few-shot GBML, we can apply gradient descent for only a limited
number of steps for fast adaptation.

Let $f(x|\theta)\in\mathbb{R}^{N}$ be the output of the classifier
parameterized by $\theta$ for a given input $x$. In this setting, $f(\cdot)$
can be considered as the logit before the softmax operation. The model
parameters $\theta$ are updated by the loss,
$L(x,y|\theta)=D(s(f(x|\theta)),y)$, where $x$ and $y$ are the input and the
corresponding label, $s(\cdot)$ is the softmax operation and $D(\cdot,\cdot)$
is some distance metric.

As shown in Fig. 1, at each inner loop, we sample a task $\tau\sim
p(\mathcal{T})$, an $N$-way classification problem consisting of two sets of
labeled data: $\mathcal{X}_{s}^{\tau}$ (support set) and
$\mathcal{X}_{t}^{\tau}$ (target set). In the inner loop, the goal is to find
the task-specific parameter $\theta^{\star}_{\tau}$ which minimizes the loss
$L$ for the support set $\mathcal{X}_{s}^{\tau}$ given the initial parameter
$\theta_{0}$ as follows:

$\text{Inner
Loop:}\quad\theta^{\star}_{\tau}=\underset{\theta}{\arg\min}\sum_{(x,y)\in\mathcal{X}_{s}^{\tau}}L(x,y|\theta;\theta_{0}).$
(1)

Normally, $\theta^{\star}_{\tau}$ is obtained by SGD (Schmidhuber, 1987; Rusu
et al., 2018a; Oh et al., 2020; Raghu et al., 2019):
$\theta^{k+1}=\theta^{k}-\alpha\nabla_{\theta}L(\theta^{k}),$ where $k$
indicates an inner loop optimization step and $\alpha$ is a learning rate. We
use $T$ steps in the inner loop.

Then, for the outer loop, we optimize the initialization parameters
$\theta_{0}$ which improve the performance of $\theta_{\tau}^{\star}$ on
$\mathcal{X}_{t}$ as follows:

$\text{Outer
Loop:}\quad\theta_{0}^{\star}=\underset{\theta_{0}}{\arg\min}\sum_{\tau}\sum_{(x,y)\in\mathcal{X}_{t}^{\tau}}L(x,y|\theta^{\star}_{\tau};\theta_{0}).$
(2)

Although GBML gained its success with this setting, this success is quite
curious: since $N$ classes constituting each task are randomly sampled and the
configuration (class) order is also random, ‘What characteristics of meta-
initialization $\theta_{0}$ make GBML work?’ remains quite uncertain.

##  4 Optimization in GBML

###  4.1 Hypothesis: Outer loop implicitly forces inner loop’s loss surface
linear

The perspective of gradient flow provides a way to understand the dynamics of
optimization. In this view, the evolution of model parameters is seen as a
discretization of continuous time evolution. Specifically, the equation
governing the model parameters is $\dot{\theta}=-\nabla_{\theta}L(\theta)$,
where $\theta$ represents the model parameter and $L(\theta)$ is the loss
function. With gradient flow, we can characterize our models into two
principal types: the Fast Model and the Slow Model. These definitions emerge
from our observation of the models’ convergence timelines. Considering the
terminal point in time as 1 and the initialization as 0, we deduce the
timestamp of a singular learning step as $\frac{1}{N}$, where N signifies the
number of steps.

Fast Model: A model with a relatively small value of N is should converge
rapidly. A salient exemplar is GBML’s inner loop, typically evolving within a
minuscule span of 1 to 5 optimization steps. Such rapid convergence
necessitates an elevated learning rate.

Slow Model: In contrast, models with a larger value of N take longer to
converge. Traditional deep learning models fall into this category. They
require many gradient steps and a lower learning rate.

This separation makes clear the differences between typical deep learning
models and GBML’s inner loop. However, we should note that GBML’s high
learning rate approach can arise some issues.

$\int_{0}^{1}\nabla L(\theta(t))\cdot\nabla
L(\theta^{k})dt\approx\int_{0}^{1}\|\nabla L(\theta^{k})\|^{2}_{2}-\alpha
t\nabla L(\theta^{k})^{T}H(\theta^{k})\nabla L(\theta^{k})dt>0,$ (3)

Gradient descent is the notion of treating the model as a (locally) linear
approximation, adjusting its parameters in a direction where the loss
decreases most rapidly. If a model adheres to Eq. 3, we can expect the
gradient descent effective in that context. (Refer to Appendix A for the
detailed proof.) However, the inner loop of GBML can be characterized as a
fast model. This implies a pronounced influence of the Hessian in the
optimization dynamics. Due to this, we cannot guarantee loss decreases in the
inner loop. In this paper, We argue that the outer loop’s primary role in GBML
is to counteract this influence of the Hessian. In more specific terms, it
aims to attenuate the impact of $H(\theta^{k})\nabla L(\theta^{k})$ on the
optimization trajectory.

Figure 3: (left) Averaged cosine similarity between the overall parameter
difference ($\theta^{\star}_{\tau}-\theta_{0}$) and each step’s gradient
$(\nabla_{\theta^{k}}L(\theta^{k}))$ in the inner loop over the training
(meta) epochs. The cosine similarity increases as training progresses,
indicating that the gradient does not vary much during an inner loop. (right)
Validation accuracy vs. training epoch. To ensure a fair comparison, we
shifted $\text{SHOT}_{p}$ for 30 epochs as it pretrains for 30 epochs with
$L_{SHOT}$ only. For this figure, we trained $\text{SHOT}_{r}$ for 300 epochs
to compare it at a glance. Our results show that SHOT converges much faster
than the baseline. Additionally, we observed a strong correlation between
cosine similarity and validation accuracy.

Fig. 3 supports our claim. As mentioned before, the optimization steps are a
discrete representation of the gradient over time. This allows us to track how
the gradient changes as time progresses. In meta-learning, the gradient’s
direction and magnitude do not change significantly as time progresses. This
is because meta-learning algorithms are designed to learn how to quickly adapt
to new tasks. By learning how to adapt to new tasks, meta-learning algorithms
can avoid the need to learn from scratch on each new task. In contrast, in
conventional deep learning, the gradient’s direction and magnitude change
significantly as time progresses. This is because conventional deep learning
algorithms are not designed to learn how to quickly adapt to new tasks.
Instead, conventional deep learning algorithms are designed to learn a general
model that can perform well on a wide variety of tasks. The difference in how
the gradient changes over time between meta-learning and conventional deep
learning is one of the key factors that explains why meta-learning algorithms
are able to learn new tasks so quickly.

It seems clear that Eq. 3 plays a crucial role in GBML, suggesting that we can
enhance its performance and dynamics by explicitly inducing the model to meet
this condition. However, the challenge remains: how can we directly reduce the
Hessian?

###  4.2 SHOT (Suppressing Hessian along the Optimization Trajectory)

In this part, we present our novel algorithm, Suppressing Hessian along the
Optimization Trajectory (SHOT), which is the main contribution of our paper.
SHOT is designed to obviate the effect of the Hessian along the inner loop
while maintaining the order of computation cost.

As previously discussed, directly reducing the Hessian is impractical.
However, we can reduce the effect of the Hessian along the optimization
trajectory, specifically by minimizing the product $H(\theta_{t})\nabla
L(\theta_{t})$, which is in line with the condition expressed in Eq. 3 and
computationally feasible. To achieve this, we define a measure that quantifies
the effect of the Hessian along the optimization trajectory.

To do so, we first create a reference model that is not influenced by the
Hessian. We obtain this reference model by increasing the number of
optimization steps. The intuition behind this is that as we increase the
number of optimization steps, the smaller the discretized step becomes, which
means that the model with many optimization steps is less influenced by the
Hessian. We then define the measure as the distance between the target model
(a model with few optimization steps, $\theta_{t}$) and the reference model (a
model with many optimization steps, $\theta_{r}$). This distance metric,
denoted as $D(\theta_{t},\theta_{r})$, quantifies the distortion along the
optimization path and evaluates the deviation from the ideal optimization
trajectory. By minimizing this measure, we can achieve our goal of reducing
the effect of the Hessian along the optimization trajectory. We will call this
method as SHOT (Suppressing Hessian along the Optimization Trajectory) which
adds the following loss term in the outer loop of a conventional GBML
algorithm:

$L_{SHOT}=D(\theta_{t}^{T},\theta_{r}^{R}).$ (4)

Here, $T$ and $R$ denote the number of steps in an inner loop used to obtain
the target and reference model respectively. In this paper, we used KL
divergence ($KL(s(f(\theta_{t}^{T}))|s(f(\theta_{r}^{R}))$) as a distance
metric $D$, but we also conducted experiments with other distance metrics such
as $L_{2}(||\theta_{t}^{T}-\theta_{r}^{R}||_{2}^{2})$ and cross-entropy
($-\sum_{c}s(f(\theta_{t}^{T}))_{c}\log s(f(\theta_{r}^{R}))_{c}$, as shown in
Table. 6. We set $T<R$ because our intention is to measure distortion caused
by the Hessian in the inner loop. As we increase the number of optimization
steps in the inner loop, the effect of distortion by the Hessian decreases.
Therefore, we could use SHOT as a pseudo-measure of distortion in the inner
loop caused by the Hessian. To match the scale of distribution between
$\theta_{t}^{T}$ and $\theta_{r}^{R}$, we set $\alpha$ in Eq. 1 as
$\alpha_{r}=\frac{T}{R}\alpha_{t}$.

Reducing the effect of Hessian along the optimization in the inner loop, as
measured by Eq. 4, can lead to faster convergence and better performance. By
minimizing the deviation from the ideal optimization trajectory, SHOT can
effectively enforce the condition expressed in Eq. 3, which is critical for
the success of GBML. As shot is the measure of distortion in the whole
trajectory of an inner loop. We use SHOT at the outer loop, where scoring of
each optimization trajectory is done.

Table 1: Comparison of the test accuracies (%) of a 4-convolutional network
with random initialization and SHOT. SHOT does not receive label information
from the outer loop. (i.e., $\text{SHOT}_{r}$ without supervision ) The
numbers in parentheses indicate the number of shots.

meta-train miniImageNet Cars meta-test miniImageNet tieredImageNet Cars Cars
CUB miniImageNet Random-init (1) 21.36 $\pm$ 0.01 21.20 $\pm$ 0.01 21.10 $\pm$
0.93 21.20 $\pm$ 0.93 21.62 $\pm$ 0.02 21.36 $\pm$ 0.01 SHOT-init (1) 24.62
$\pm$ 0.04 23.80 $\pm$ 0.06 24.63 $\pm$ 0.04 25.84 $\pm$ 0.06 25.86 $\pm$ 0.03
26.48 $\pm$ 0.04 Random-init (5) 21.29 $\pm$ 0.09 21.72 $\pm$ 0.19 21.58 $\pm$
0.37 21.58 $\pm$ 0.04 22.31 $\pm$ 0.06 21.28 $\pm$ 0.09 SHOT-init (5) 27.32
$\pm$ 0.59 31.28 $\pm$ 1.68 29.10 $\pm$ 1.81 28.28 $\pm$ 0.79 30.34 $\pm$ 2.69
26.61 $\pm$ 1.05

Table 1 highlights the importance of distortion in GBML. We compared the
performance of a randomly-initialized model with that of a lower-distortion-
initialized model, where we only reduced the distortion introduced in Eq. 4 in
the outer loop. As our SHOT loss is self-supervisory that can be applied
without target class labels, no target labels are given here. The results
demonstrate that the lower-distortion-initialized model outperforms the
randomly-initialized model, indicating that distortion is a crucial factor in
GBML. This also suggests that we can generally improve GBML by incorporating
our loss term, SHOT.

###  4.3 SHOT can be used as a regularizer or meta-meta-initializer in the
outer loop.

One potential approach to incorporating SHOT is to use it as a regularizer in
GBML. This approach involves adding SHOT to the loss function of the outer
loop during training, thereby leveraging it as a regularizer reducing the
effect of the Hessian along the optimization trajectory. Alternatively, SHOT
can be utilized as a meta-meta-initializer as is used to obtain Table 1. In
this approach, the entire model is pre-trained with SHOT, i.e., The outer loss
is substituted with $L_{SHOT}$ without utilizing any labeling information.
Then the model is trained with the original loss function. However, there is a
vulnerability to homogeneous corruption, where the loss goes to zero only if
$\theta_{t}$ and $\theta_{r}$ produce the same output of all inputs. As this
phenomenon is trivial in the self-supervised learning area, we drew
inspiration from BYOL Grill et al. (2020). To address this issue, we added a
projector to the reference model, while the target remained the same. We also
augmented $\mathcal{X}_{t}^{\tau}$. After pre-training with SHOT, we started
training by loading the best model with the highest validation accuracy. To
differentiate between SHOT utilized as a regularizer and SHOT utilized as a
meta-meta-initializer, we denote them as $\text{SHOT}_{r}$ and
$\text{SHOT}_{p}$, respectively.

###  4.4 SHOT does not require any computation burden

Several meta-learning algorithms have already assumed linearity in the inner
loop. For example, FoMAML and Reptile Finn et al. (2017); Nichol et al. (2018)
proposed Hessian-Free GBML algorithms that both assumed linearity in the inner
loop. Additionally, ANIL Raghu et al. (2019) proposed a feature reuse concept
that allows for translation to zero distortion in the inner loop. However,
none of these studies explicitly argued that GBML implicitly suppresses the
Hessian in the inner loop as an inductive prior, nor did they propose an
algorithm that controls the Hessian. Our algorithm is similar to the ‘SAM’
optimizer Foret et al. (2020) or natural gradient descent Bonnabel (2013).
Like SAM and natural gradient descent, our algorithm searches for a good
initialization point based on its curvature. However, SHOT finds it without
much increasing the computational cost compared to normal gradient descent,
which is a notable advantage over SAM or natural gradient descent. Although
SAM uses a first-order computation term, it searches for flat minima by
iteratively adding noise to each point. Natural gradient descent exploits the
Fisher metric, which is a second-order computation cost.

On the other hand, we can make SHOT not give any computation burden. When we
set the reference model $\theta_{r}$ to exploit the same number of
optimization steps in the inner loop as the baseline model and the target
model $\theta_{t}$ exploits one optimization step in the inner loop, because
both models share the same initialization point $\theta^{0}$, i.e.,
$\theta_{t}^{0}=\theta_{r}^{0}$. $\theta_{t}$ can be obtained by calculating
the first optimization step of $\theta_{r}$. Therefore, SHOT does not require
any additional backpropagation. What we require is only one more inference at
the endpoint $\theta_{t}^{T}$ for ${\text{SHOT}_{r}}$. Also, even if we set
$\theta_{r}$ to have more optimization steps than the baseline model and set
$\theta_{t}$ to have the same number of optimization steps as the baseline
model, the computation cost at test time is unchanged as it is only used in
the training phase.

##  5 Experiments

Datasets  Our method is evaluated on miniImageNet (Vinyals et al., 2016),
tieredImageNet Ren et al. (2018), Cars Krause et al. (2013), and CUB Welinder
et al. (2010) datasets. miniImageNet and tieredImageNet are subsets of
ImageNet Russakovsky et al. (2015) with 100 classes each and can be considered
as general datasets. Cars and CUB are widely used benchmarks for fine-grained
image classification tasks as they contain similar objects with subtle
differences. By conducting experiments on these datasets, we can evaluate the
performance of our method across a range of tasks.

Architectures  We use 4-Conv from Finn et al. (2017) and ResNet-12 from
Oreshkin et al. (2018) as our backbone architecture. We trained 4-Conv for 30k
iterations and ResNet-12 for 10k iterations.

Experimental Setup  To ensure a fair comparison, we fine-tuned MAML Finn et
al. (2017) first, and our baselines outperformed the original paper. We used
Adam Kingma and Ba (2014) as our optimizer and fixed the outer loop learning
rate to $10^{-3}$ and inner loop learning rate to 0.5 for 4-Conv, and
$10^{-5}$ and 0.01 for ResNet-12. We restricted the total number of steps in
the inner loop to 3, as argued in Finn et al. (2017). We set the number of
steps to 1 for the target model $\theta_{t}$ and 3 for the reference model
$\theta_{r}$ in our SHOT. With this setting, we can exploit SHOT without any
additional computation burden as discussed in Sec 4.4. For $\text{SHOT}_{r}$,
we fixed $\lambda$ to 0.1 when we set the original model as $\theta_{r}$. For
one-step algorithms, we set the number of optimization steps to 2 for
$\theta_{r}$ and 1 for $\theta_{t}$. So actual optimization step in inference
time is 1, which is the same as the baseline. We set $\lambda$ for $10^{-6}$
in that case. For the Hessian-Free algorithm FoMAML, we fine-tuned $\lambda$
to meet the best performance. We set the learning rate of the outer and inner
loops to $10^{-4}$ and 0.01, respectively, for 4-Conv. To maintain
conciseness, we only report the results of $\text{SHOT}_{r}$ except for the
comparison with MAML baseline. We conducted our experiments on a single A100
GPU. We implemented our method using torchmeta Deleu et al. (2019) library
except for the implementation of LEO Rusu et al. (2018a), where we used the
official implementation and embeddings.

Table 2: Test accuracy % of 4-conv network on benchmark data sets. The values
in parentheses indicate the number of shots. The best accuracy among different
methods is bold-faced.

Domain General (Coarse-grained) Specific (Fine-grained) Dataset miniImageNet
tieredImageNet Cars Cub MAML (1) 47.88 $\pm$ 0.55 46.93 $\pm$ 1.07 47.78 $\pm$
0.99 57.04 $\pm$ 1.42 MAML + $\text{SHOT}_{r}$ (1) 47.97 $\pm$ 0.71 47.53
$\pm$ 0.59 50.44 $\pm$ 0.62 57.55 $\pm$ 0.64 MAML + $\text{SHOT}_{p}$ (1)
48.11 $\pm$ 0.26 47.53 $\pm$ 0.68 49.08 $\pm$ 0.88 58.23 $\pm$ 0.39 MAML (5)
64.81 $\pm$ 1.63 66.12 $\pm$ 1.10 62.24 $\pm$ 2.01 72.48 $\pm$ 0.86 MAML +
$\text{SHOT}_{r}$ (5) 66.86 $\pm$ 0.58 69.08 $\pm$ 0.31 64.20 $\pm$ 1.58 73.38
$\pm$ 0.32 MAML + $\text{SHOT}_{p}$ (5) 66.35 $\pm$ 0.27 68.94 $\pm$ 0.87
64.84 $\pm$ 2.87 73.27 $\pm$ 0.40

Table 3: Test accuracy % of 4-conv network on cross-domain adaptation. The
values in parentheses indicate the number of shots. The best accuracy among
different methods is bold-faced.

adaptation General to General General to Specific Specific to General Specific
to Specific meta-train tieredImageNet miniImageNet miniImageNet miniImageNet
Cars Cars CUB Cars meta-test miniImageNet tieredImageNet Cars CUB miniImageNet
tieredImageNet Cars CUB MAML (1) 47.52 $\pm$ 1.66 51.84 $\pm$ 0.24 34.41 $\pm$
0.47 40.91 $\pm$ 0.57 28.67 $\pm$ 1.17 30.79 $\pm$ 1.17 32.74 $\pm$ 1.12 30.95
$\pm$ 1.41 MAML + $\text{SHOT}_{r}$ (1) 48.20 $\pm$ 0.64 51.68 $\pm$ 0.68
34.03 $\pm$ 1.11 41.58 $\pm$ 0.56 30.19 $\pm$ 0.50 31.96 $\pm$ 0.50 33.11
$\pm$ 0.41 31.21 $\pm$ 0.64 MAML + $\text{SHOT}_{p}$ (1) 48.45 $\pm$ 0.13
51.97 $\pm$ 0.25 34.84 $\pm$ 0.16 41.65 $\pm$ 0.70 29.02 $\pm$ 0.73 31.26
$\pm$ 0.54 33.69 $\pm$ 0.45 31.65 $\pm$ 0.39 MAML (5) 66.86 $\pm$ 1.58 67.96
$\pm$ 1.22 46.57 $\pm$ 0.53 56.32 $\pm$ 1.17 37.23 $\pm$ 1.95 41.00 $\pm$ 1.86
44.02 $\pm$ 2.29 41.84 $\pm$ 1.25 MAML + $\text{SHOT}_{r}$ (5) 70.70 $\pm$
0.29 69.48 $\pm$ 0.17 48.42 $\pm$ 0.72 58.40 $\pm$ 0.48 40.79 $\pm$ 0.93 42.73
$\pm$ 0.32 44.47 $\pm$ 0.25 43.46 $\pm$ 0.89 MAML + $\text{SHOT}_{p}$ (5)
69.88 $\pm$ 0.53 69.83 $\pm$ 0.34 47.99 $\pm$ 2.10 58.30 $\pm$ 0.44 37.96
$\pm$ 1.14 41.94 $\pm$ 0.23 43.83 $\pm$ 1.06 40.64 $\pm$ 1.25

Table 4: Test accuracy % of ResNet-12. The values in parentheses indicate the number of shots. The better accuracy between the baseline and SHOT is bold-faced. meta-train | miniImageNet
---|---
meta-test | miniImageNet | tieredImageNet | Cars
MAML (1) |  49.47 $\pm$ 0.43 |  54.88 $\pm$ 1.07 |  32.74 $\pm$ 0.48
MAML + $\text{SHOT}_{r}$ (1) |  51.11 $\pm$ 0.32 |  55.14 $\pm$ 0.63 |  32.66 $\pm$ 0.27
MAML (5) |  69.12 $\pm$ 0.08 |  71.60 $\pm$ 0.08 |  51.99 $\pm$ 0.25
MAML + $\text{SHOT}_{r}$ (5) |  69.74 $\pm$ 0.16 |  71.71 $\pm$ 0.39 |  52.16 $\pm$ 0.28

###  5.1 Characteristics of SHOT

SHOT Outperforms the Baselines on Benchmarks and Cross-Benchmark Settings  Our
method, SHOT, outperforms the corresponding baseline on all benchmarks, as
shown in Table 2. Both $\text{SHOT}_{r}$ and $\text{SHOT}_{p}$ outperform the
baseline, regardless of whether the benchmark is general or fine-grained.
Additionally, SHOT generalizes well to cross-benchmark settings, as
demonstrated in Table 3, where it improves the performance of the baseline.
Also, we can see that this result is independent of the backbone architecture,
as shown in Table 4.

SHOT works with Hessian-Free and One-Step algorithms  As many GBML algorithms
are based on Hessian-Free algorithms which only exploit gradients and One-Step
algorithms which only exploit one step in the inner loop Finn et al. (2017);
Nichol et al. (2018); Oh et al. (2020); Raghu et al. (2019), it is important
to show that our method also works with these algorithms. As shown in Table 5,
our method improves the performance of the baseline on both Hessian-Free and
One-Step algorithms. For One-Step algorithms, we set $\theta_{t}$ as the
original model and $\theta_{r}$ as the 2-step optimized model. This is
practically useful since many algorithms exploit only one step in the inner
loop. Notably, our method is applied only in the training phase, so it has the
same inference cost as the baseline. Also, we only used first-order
optimization for the meta-optimizer to ensure a fair comparison when adding
SHOT to FoMAML. For ANIL Raghu et al. (2019) and BOIL Oh et al. (2020), we
conducted experiments with the same settings from each baseline paper, and we
bring the results from Oh et al. (2020) for comparison.

Table 5: Test accuracy (%) of 4-conv network on benchmark data sets. The
values in parentheses indicate the number of shots. The better accuracy
between the baseline and SHOT is bold-faced. Baselines in the table
corresponds to either one-step algorithm (which uses one step optimization in
the inner loop) or Hessian-Free algorithm (an algorithm which only uses
gradients in the outer loop).

Domain General(Coarse-grained) Specific (Fine-grained) Dataset miniImageNet
tieredImageNet Cars Cub FoMAML (1) 47.70 $\pm$ 0.56 47.72 $\pm$ 0.76 49.08
$\pm$ 1.46 60.18 $\pm$ 0.45 FoMAML + $\text{SHOT}_{r}$ (1) 48.49 $\pm$ 0.29
47.91 $\pm$ 0.59 48.59 $\pm$ 0.72 60.54 $\pm$ 0.29 ANIL (1) 47.82 $\pm$ 0.20
49.35 $\pm$ 0.26 46.81 $\pm$ 0.24 57.03 $\pm$ 0.41 ANIL + $\text{SHOT}_{r}$
(1) 48.02 $\pm$ 0.38 48.29 $\pm$ 0.93 47.82 $\pm$ 1.48 57.31 $\pm$ 0.89 BOIL
(1) 49.61 $\pm$ 0.16 48.58 $\pm$ 0.27 56.82 $\pm$ 0.21 61.60 $\pm$ 0.57 BOIL +
$\text{SHOT}_{r}$ (1) 50.36 $\pm$ 0.42 50.36 $\pm$ 0.49 56.22 $\pm$ 1.96 62.43
$\pm$ 0.30 FoMAML (5) 64.49 $\pm$ 0.46 65.25 $\pm$ 0.64 67.99 $\pm$ 1.20 73.18
$\pm$ 0.18 FoMAML + $\text{SHOT}_{r}$ (5) 64.58 $\pm$ 0.06 65.92 $\pm$ 0.22
68.02 $\pm$ 1.00 73.89 $\pm$ 0.18 ANIL (5) 63.04 $\pm$ 0.42 65.82 $\pm$ 0.12
61.95 $\pm$ 0.38 70.93 $\pm$ 0.28 ANIL + $\text{SHOT}_{r}$ (5) 64.40 $\pm$
0.75 66.19 $\pm$ 0.13 61.00 $\pm$ 0.90 71.71 $\pm$ 1.34 BOIL (5) 66.45 $\pm$
0.37 69.37 $\pm$ 0.12 75.18 $\pm$ 0.21 75.96 $\pm$ 0.17 BOIL +
$\text{SHOT}_{r}$ (5) 66.34 $\pm$ 0.21 69.20 $\pm$ 0.25 75.25 $\pm$ 0.14 76.77
$\pm$ 0.51

SHOT enables faster convergence. Our method, SHOT, achieves faster convergence
compared to the baseline, as demonstrated in Figure 3, where we plot the
validation accuracy against the training epoch. To understand why this occurs,
we plotted averaged cosine similarity between the overall parameter difference
with each inner gradient step in the inner loop, as shown in Figure 3. We
found that overall similarity increases in both cases of with and without
SHOT, and the similarity is highly correlated with the validation accuracy.
SHOT enforces high similarity from the beginning, resulting in faster
convergence with better generalization ability.

Applying other distance metrics. We used KL divergence as a distance metric in
SHOT, but it is possible to use other metrics as well. For example, we could
use cross-entropy as a distance metric from a probability perspective, or we
could use L2 distance between parameters as a distance metric. We conducted
experiments with different distance metrics, and for L2 distance we divided
the L2 distance by the square root of the number of parameters to account for
the initialization scheme Kumar (2017). For the cross-entropy loss, we used
the same $\lambda$ (balancing term for Eq. 4) settings as for KL divergence.
The results are shown in Table 6. We observed that SHOT outperformed the
baseline for every distance metric. This supports our hypothesis that we can
improve performance and dynamics by suppressing the effect of the Hessian
along the optimization trajectory.

Table 6: Test accuracy % of 4-conv using SHOT with various distance metrics. The values in parentheses indicate the number of shots. The better accuracy between the baseline and SHOT is bold-faced. meta-train | miniImageNet
---|---
meta-test | miniImageNet | tieredImageNet | Cars
Baseline (1) |  47.88 $\pm$ 0.55 |  51.84 $\pm$ 0.24 |  34.41 $\pm$ 0.47
KL-Divergence (1) |  47.97 $\pm$ 0.71 |  51.68 $\pm$ 0.68 |  34.03 $\pm$ 1.11
Cross-Entropy (1) |  49.07 $\pm$ 0.14 |  53.83 $\pm$ 0.14 |  34.54 $\pm$ 0.31
L2 Distance (1) |  50.89 $\pm$ 0.06 |  54.63 $\pm$ 0.25 |  35.75 $\pm$ 0.29
Baseline (5) |  64.81 $\pm$ 1.63 |  67.96 $\pm$ 1.22 |  46.57 $\pm$ 0.53
KL-Divergence (5) |  66.86 $\pm$ 0.58 |  71.77 $\pm$ 0.42 |  48.42 $\pm$ 0.72
Cross-Entropy (5) |  66.34 $\pm$ 0.35 |  69.47 $\pm$ 0.35 |  48.86 $\pm$ 0.71
L2 Distance (5) |  66.85 $\pm$ 0.44 |  69.96 $\pm$ 0.27 |  48.40 $\pm$ 0.12

##  6 Conclusion

In this paper, we claim that GBML implicitly suppresses the Hessian along the
optimization path. We also demonstrate that we can improve the performance of
GBML by explicitly enforcing this implicit prior. To achieve this, we propose
a novel method, SHOT (Suppressing the Hessian along Optimization Trajectory).
Although it deals with curvature of the loss surface, it has the same
computational order compared to other GBML methods. More specifically, it adds
a negligible computational overhead to the baseline and costs the same at
inference time. We show that SHOT outperforms the baseline and is algorithm-
and architecture-independent. We believe that our method can be applied to
other GBML methods in a plug-and-play manner.

##  7 Acknowledgements

This work was supported by NRF grant (2021R1A2C3006659) and IITP grants
(2021-0-01343, 2022-0-00953), all of which were funded by Korean Government
(MSIT).

## References

  * Baik et al. [2020] S. Baik, M. Choi, J. Choi, H. Kim, and K. M. Lee.  Meta-learning with adaptive hyperparameters.  _Advances in Neural Information Processing Systems_ , 33:20755–20765, 2020.
  * Bernacchia [2021] A. Bernacchia.  Meta-learning with negative learning rates.  _arXiv preprint arXiv:2102.00940_ , 2021.
  * Bonnabel [2013] S. Bonnabel.  Stochastic gradient descent on riemannian manifolds.  _IEEE Transactions on Automatic Control_ , 58(9):2217–2229, 2013.
  * Bromley et al. [1993] J. Bromley, I. Guyon, Y. LeCun, E. Säckinger, and R. Shah.  Signature verification using a" siamese" time delay neural network.  _Advances in neural information processing systems_ , 6, 1993.
  * Chen et al. [2020] T. Chen, S. Kornblith, M. Norouzi, and G. Hinton.  A simple framework for contrastive learning of visual representations.  In _International conference on machine learning_ , pages 1597–1607. PMLR, 2020.
  * Chen and He [2021] X. Chen and K. He.  Exploring simple siamese representation learning.  In _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , pages 15750–15758, 2021.
  * Deleu et al. [2019] T. Deleu, T. Würfl, M. Samiei, J. P. Cohen, and Y. Bengio.  Torchmeta: A Meta-Learning library for PyTorch, 2019.  URL https://arxiv.org/abs/1909.06576.  Available at: https://github.com/tristandeleu/pytorch-meta.
  * Dosovitskiy et al. [2020] A. Dosovitskiy, L. Beyer, A. Kolesnikov, D. Weissenborn, X. Zhai, T. Unterthiner, M. Dehghani, M. Minderer, G. Heigold, S. Gelly, J. Uszkoreit, and N. Houlsby.  An image is worth 16x16 words: Transformers for image recognition at scale.  _CoRR_ , abs/2010.11929, 2020.  URL https://arxiv.org/abs/2010.11929.
  * Fallah et al. [2020] A. Fallah, A. Mokhtari, and A. Ozdaglar.  On the convergence theory of gradient-based model-agnostic meta-learning algorithms.  In S. Chiappa and R. Calandra, editors, _Proceedings of the Twenty Third International Conference on Artificial Intelligence and Statistics_ , volume 108 of _Proceedings of Machine Learning Research_ , pages 1082–1092. PMLR, 26–28 Aug 2020.  URL https://proceedings.mlr.press/v108/fallah20a.html.
  * Finn et al. [2017] C. Finn, P. Abbeel, and S. Levine.  Model-agnostic meta-learning for fast adaptation of deep networks.  In _International conference on machine learning_ , pages 1126–1135. PMLR, 2017.
  * Flennerhag et al. [2019] S. Flennerhag, A. A. Rusu, R. Pascanu, F. Visin, H. Yin, and R. Hadsell.  Meta-learning with warped gradient descent.  _arXiv preprint arXiv:1909.00025_ , 2019.
  * Foret et al. [2020] P. Foret, A. Kleiner, H. Mobahi, and B. Neyshabur.  Sharpness-aware minimization for efficiently improving generalization.  _CoRR_ , abs/2010.01412, 2020.  URL https://arxiv.org/abs/2010.01412.
  * Grill et al. [2020] J.-B. Grill, F. Strub, F. Altché, C. Tallec, P. Richemond, E. Buchatskaya, C. Doersch, B. Avila Pires, Z. Guo, M. Gheshlaghi Azar, et al.  Bootstrap your own latent-a new approach to self-supervised learning.  _Advances in neural information processing systems_ , 33:21271–21284, 2020.
  * He et al. [2018] T. He, Z. Zhang, H. Zhang, Z. Zhang, J. Xie, and M. Li.  Bag of tricks for image classification with convolutional neural networks, 2018.  URL https://arxiv.org/abs/1812.01187.
  * Hiller et al. [2022] M. Hiller, M. Harandi, and T. Drummond.  On enforcing better conditioned meta-learning for rapid few-shot adaptation.  In _NeurIPS_ , 2022.
  * Hinton et al. [2015] G. Hinton, O. Vinyals, and J. Dean.  Distilling the knowledge in a neural network.  _arXiv preprint arXiv:1503.02531_ , 2015.
  * Kingma and Ba [2014] D. P. Kingma and J. Ba.  Adam: A method for stochastic optimization.  _arXiv preprint arXiv:1412.6980_ , 2014.
  * Krause et al. [2013] J. Krause, M. Stark, J. Deng, and L. Fei-Fei.  3d object representations for fine-grained categorization.  _In Proceedings of the IEEE international conference on computer vision workshops_ , 2013.
  * Kumar [2017] S. K. Kumar.  On weight initialization in deep neural networks.  _arXiv preprint arXiv:1704.08863_ , 2017.
  * LeCun et al. [2015] Y. LeCun, Y. Bengio, and G. Hinton.  Deep learning.  _nature_ , 521(7553):436–444, 2015.
  * Nichol et al. [2018] A. Nichol, J. Achiam, and J. Schulman.  On first-order meta-learning algorithms.  _arXiv preprint arXiv:1803.02999_ , 2018.
  * Oh et al. [2020] J. Oh, H. Yoo, C. Kim, and S.-Y. Yun.  Boil: Towards representation change for few-shot learning.  _arXiv preprint arXiv:2008.08882_ , 2020.
  * Oreshkin et al. [2018] B. N. Oreshkin, P. R. López, and A. Lacoste.  TADAM: task dependent adaptive metric for improved few-shot learning.  _CoRR_ , abs/1805.10123, 2018.  URL http://arxiv.org/abs/1805.10123.
  * Park and Oliva [2019] E. Park and J. B. Oliva.  Meta-curvature.  In _NeurIPS_ , 2019.  Published 9 February 2019.
  * Radford et al. [2021] A. Radford, J. W. Kim, C. Hallacy, A. Ramesh, G. Goh, S. Agarwal, G. Sastry, A. Askell, P. Mishkin, J. Clark, G. Krueger, and I. Sutskever.  Learning transferable visual models from natural language supervision.  _CoRR_ , abs/2103.00020, 2021.  URL https://arxiv.org/abs/2103.00020.
  * Raghu et al. [2019] A. Raghu, M. Raghu, S. Bengio, and O. Vinyals.  Rapid learning or feature reuse? towards understanding the effectiveness of maml.  _arXiv preprint arXiv:1909.09157_ , 2019.
  * Rajeswaran et al. [2019] A. Rajeswaran, C. Finn, S. M. Kakade, and S. Levine.  Meta-learning with implicit gradients.  _Advances in neural information processing systems_ , 32, 2019.
  * Ravi and Larochelle [2017] S. Ravi and H. Larochelle.  Optimization as a model for few-shot learning.  In _International Conference on Learning Representations_ , 2017.  URL https://openreview.net/forum?id=rJY0-Kcll.
  * Ren et al. [2018] M. Ren, E. Triantafillou, S. Ravi, J. Snell, K. Swersky, J. B. Tenenbaum, H. Larochelle, and R. S. Zemel.  Meta-learning for semi-supervised few-shot classification.  _In the Sixth International Conference on Learning Representations_ , 2018.
  * Russakovsky et al. [2015] O. Russakovsky, J. Deng, H. Su, J. Krause, S. Satheesh, S. Ma, Z. Huang, A. Karpathy, A. Khosla, M. Bernstein, et al.  Imagenet large scale visual recognition challenge.  _International journal of computer vision_ , 115:211–252, 2015.
  * Rusu et al. [2018a] A. A. Rusu, D. Rao, J. Sygnowski, O. Vinyals, R. Pascanu, S. Osindero, and R. Hadsell.  Meta-learning with latent embedding optimization.  _CoRR_ , abs/1807.05960, 2018a.  URL http://arxiv.org/abs/1807.05960.
  * Rusu et al. [2018b] A. A. Rusu, D. Rao, J. Sygnowski, O. Vinyals, R. Pascanu, S. Osindero, and R. Hadsell.  Meta-learning with latent embedding optimization.  _arXiv preprint arXiv:1807.05960_ , 2018b.
  * Schmidhuber [1987] J. Schmidhuber.  _Evolutionary principles in self-referential learning, or on learning how to learn: the meta-meta-… hook_.  PhD thesis, Technische Universität München, 1987.
  * Simon et al. [2020] C. Simon, P. Koniusz, R. Nock, and M. Harandi.  On modulating the gradient for meta-learning.  In _ECCV_ , 2020.
  * Snell et al. [2017] J. Snell, K. Swersky, and R. Zemel.  Prototypical networks for few-shot learning.  _Advances in neural information processing systems_ , 30, 2017.
  * Sung et al. [2018] F. Sung, Y. Yang, L. Zhang, T. Xiang, P. H. Torr, and T. M. Hospedales.  Learning to compare: Relation network for few-shot learning.  In _Proceedings of the IEEE conference on computer vision and pattern recognition_ , pages 1199–1208, 2018.
  * Vinyals et al. [2016] O. Vinyals, C. Blundell, T. Lillicrap, and D. Wierstra.  Matching networks for one shot learning.  _Advances in neural information processing systems_ , 29, 2016.
  * Wang et al. [2020] H. Wang, R. Sun, and B. Li.  Global convergence and induced kernels of gradient-based meta-learning with neural nets.  _CoRR_ , abs/2006.14606, 2020.  URL https://arxiv.org/abs/2006.14606.
  * Welinder et al. [2010] P. Welinder, S. Branson, T. Mita, C. Wah, F. Schroff, S. Belongie, and P. Perona.  Caltech-ucsd birds 200.  _Technical Report CNS-TR-2010-001, California Institute of Technology_ , 2010.

##  Appendix A Proof of (3)

Due to fundamental theorem of calculus,

$L(\theta^{k+1})-L(\theta^{k})=\int_{C}\nabla L(\theta)\cdot
d\theta=\int_{0}^{1}\nabla L(\theta(t))\cdot v(t)dt,$ (5)

where $C$ is a trajectory whose start and end points are $\theta^{k}$ and
$\theta^{k+1}$. In GD setting, because $\theta^{k+1}=\theta^{k}-\alpha\nabla
L(\theta^{k})$ for some learning rate $\alpha>0$, we can think of the straight
line trajectory joining $\theta^{k}$ and $\theta^{k+1}$. In this case, the
velocity vector becomes $v(t)=-\alpha\nabla L(\theta^{k})$ and

$L(\theta^{k+1})-L(\theta^{k})=-\alpha\int_{0}^{1}\nabla
L(\theta(t))\cdot\nabla L(\theta^{k})dt$ (6)

where $\theta(0)=\theta^{k}$, $\theta(1)=\theta^{k+1}$ and
$\theta(t)=(1-t)\theta(0)+t\theta(1)$.

By Taylor series expansion it becomes

$\nabla L(\theta(t))\approx\nabla
L(\theta^{k})+H(\theta^{k})(\theta(t)-\theta^{k})=(I-\alpha
tH(\theta^{k}))\nabla L(\theta^{k})$ (7)

and combining Eq.(6) and Eq.(7) proves Eq.(3).

##  Appendix B SHOT is robust to hyperparameter settings

Table 7: Test accuracy % of 4-conv using SHOT with different $\lambda$ values. The values in parentheses indicate the number of shots. The better accuracy between the baseline and SHOT is bold-faced. size of $\lambda$ | miniImageNet (5)
---|---
0 (baseline, MAML) | 64.81 $\pm$ 1.63
0.1 ($\text{SHOT}_{r}$) |  66.86 $\pm$ 0.58
$10^{-2}$ | 65.97 $\pm$ 0.93
$10^{-3}$ | 66.85 $\pm$ 0.47
$10^{-4}$ | 66.08 $\pm$ 0.36
$10^{-5}$ | 66.81 $\pm$ 0.28

Table 7 shows the performance of SHOT with various hyperparameter settings. We
conducted experiments with different learning rates of $\lambda$, and in all
cases, SHOT improved the performance of the baseline. This suggests that SHOT
is robust to hyperparameter settings.

##  Appendix C Another viewpoint of ANIL and BOIL

In this section, we reconcile the opposite opinions of feature reuse versus
feature adaptation in gradient-based meta-learning (GBML) with our hypothesis
that reducing the impact of the Hessian in the inner loop can improve
performance.

ANIL, proposed in Raghu et al. [2019], argues that feature reuse is key in the
inner loop, where the feature remains invariant while only the decision
boundary is adapted. On the other hand, BOIL, proposed in Oh et al. [2020],
argues that feature adaptation is key in the inner loop, where the feature is
adapted while the decision boundary remains invariant.

To test their arguments, ANIL and BOIL proposed two algorithms. ANIL freezes
the encoder and only updates the head in the inner loop, while BOIL freezes
the head and only updates the encoder in the inner loop. The problem is that
both algorithm shown good performance thereby both arguments look persuasive
despite they argue exactly in the opposite ways.

Our hypothesis, which suggests that the outer loop implicitly suppresses the
Hessian along the optimization trajectory, can reconcile the arguments of both
ANIL and BOIL. This is because our hypothesis implies that the model acts
linearly in the inner loop. ANIL and BOIL can be interpreted as algorithms
that enforce linearity in the inner loop by restricting parameters and
reducing the number of non-linear components between layers.

ANIL freezes the encoder and only updates the head in the inner loop, reducing
the number of non-linear components in the inner loop. This enforces linearity
in the inner loop, as the only non-linearity is the loss function. ANIL
achieves better performance than MAML in 1-step optimization, as it is more
powerful at 1-step optimization, which views the model as linear.

BOIL freezes the head and only updates the encoder in the inner loop, reducing
the number of non-linear components in the inner loop. By applying BOIL, the
gradient norm is predominant in the last layer of the encoder, making it a
variant of ANIL that updates only the penultimate layer. This layer has much
stronger performance, as it can change the feature while maintaining the
linearity of the model. Table 14 of Oh et al. [2020] shows a boosted
performance when all but the penultimate layer is not frozen. By explicitly
enforcing linearity in the inner loop, BOIL achieves improved performance.

##  Appendix D GBML is a variant of Prototype Vector method

In this section, we provide a novel viewpoint of GBML, that GBML (Gradient-
Based Meta Learning) is a varient of MBML (Metric-Based Meta Learning). This
viewpoint relies on the linearity assumption. i.e., the effect of the Hessian
along the optimization trajectory is zero, thereby the model act as linear in
the inner loop.

Suppose there exists a meta-learning model that satisfies the linearity
assumption in the inner loop, then classifying a new classification task with
a task-specific function $f(\cdot|\theta^{\star})$ after an inner loop is
equivalent to creating a prototype vector for each class on a specific feature
map and classifying the input as the class of the most similar prototype
vector.

The proof starts by defining the prototype vector at first.

Prototype Vector We define a prototype vector $V_{c}$ for class $c$ in an
$N$-way $K$-shot classification task formally as

$V_{c}=\sum_{i=1}^{N}\sum_{j=1}^{K}\beta_{ij}\varphi_{c}(X_{ij}),\quad
c\in\\{1,\cdots,N\\},$ (8)

where $X_{ij}$ is the $j$-th input sample for the $i$-th class,
$\varphi_{c}(\cdot)\in\mathcal{H}$ is a class-specific feature map and
$\beta_{ij}$ indicates the importance of $X_{ij}$ for constituting the
prototype vector $V_{c}$.

In other words, there exists a feature map $\varphi_{c}$ for each class $c$,
and the support set is mapped to the corresponding feature map and then
weighted-averaged to constitute the prototype vector of the corresponding
class. At inference time, the classification of a given query $X$ is done by
taking the class of the most similar prototype vector as follows:

$\hat{c}=\underset{c}{\arg\max}\langle{V_{c},\varphi(X)}\rangle.$ (9)

Here, $\varphi:\mathcal{X}\rightarrow H$ is a non-class-specific mapping. We
can also rewrite the prototype vector using $\varphi$ and by defining a
projection $P_{c}:\mathcal{X}\rightarrow\mathcal{X}$ as

$\displaystyle P_{c}(X)=\begin{dcases*}X&if $y(X)=c$,\\\
\nu\in\mathcal{N}(\varphi),&if $y(X)\neq c$\end{dcases*}$ (10)

where $y(X)$ is the ground truth class of $X$ and $\mathcal{N}$ is the null
space of $\varphi$ i.e., $\varphi(\nu)=0$.

Then by defining $\varphi_{c}\triangleq\varphi\circ P_{c}$ and
$\beta_{ij}\triangleq\frac{1}{K}$, it becomes

$V_{c}=\frac{1}{K}\sum_{j=1}^{K}\varphi(X_{cj}).$ (11)

SGD in the inner loop If GBML satisfies the hypothesis of linearity in the
inner loop, $f$ is locally linear in $\theta$ in an inner loop. More
specifically, there exists an equivalent feature map
$\varphi_{c}:\mathcal{X}\rightarrow\mathcal{H}$ which satisfies
$f_{c}(\cdot|\theta_{c})=\langle\theta_{c},\varphi_{c}(\cdot)\rangle$ for
every $x\in\mathcal{X}$ where
$f(\cdot|\theta)=[f_{1}(\cdot|\theta_{1}),\cdots,f_{N}(\cdot|\theta_{N})]^{T}$.
With the loss function $L(x,y|\theta)=D(s(f(x|\theta)),y)$ for some distance
measure $D$ such as cross entropy, we can formulate the inner loop of $N$-way
$K$-shot meta learning by SGD as

$\theta_{c}^{k+1}=\theta^{k}_{c}-\alpha\sum_{i=1}^{N}\sum_{j=1}^{K}\frac{\partial
L(X_{ij},y(X_{ij})|\theta)}{\partial\theta_{c}}=\theta^{k}_{c}-\alpha\sum_{i=1}^{N}\sum_{j=1}^{K}\frac{\partial
D}{\partial f_{c}}\varphi_{c}(X_{ij}),$ (12)

since all samples in the support set are inputted in a batch of an inner loop.

Because the model is linear in the inner loop, the batch gradient does not
change. Let $\beta_{ij}=-\frac{\partial D}{\partial
f_{c}}\rvert_{\theta_{c}^{0},X_{ij}}$. Then after $t$ steps, by (8), the model
becomes

$\theta^{t}_{c}=\theta^{0}_{c}+\alpha
t\sum_{i=1}^{N}\sum_{j=1}^{K}\beta_{ij}\varphi_{c}(X_{ij})=\theta^{0}_{c}+\alpha
tV_{c}.$ (13)

At the initialization step of an inner loop, there is no information about the
class, even the configuration order of the class, because the task is randomly
sampled. If so, the problem is solved in the inner loop. For example, if a
class Dog is allocated to a specific index such as Class 3. There is no
guarantee that it will have the identical index the next time the class Dog
comes in. Thus, at a meta-initialization point $\theta^{0}$, the scores for
different classes would not be much different, i.e.,
$f_{i}(x|\theta^{0}_{i})\simeq f_{j}(x|\theta^{0}_{j})$ for
$i,j\in[1,\cdots\,N]$.

Considering the goal of classification is achieved through relative values
between $f_{i}(X)$’s, the value at the initialization point does not need to
be considered significantly. Therefore

$\underset{c}{\arg\max}f_{c}(X)=\underset{c}{\arg\max}\langle\theta^{t}_{c},\varphi(X)\rangle=\underset{c}{\arg\max}\langle\theta_{c}^{0}+\alpha
tV_{c},\varphi(X)\rangle\sim\underset{c}{\arg\max}\langle
V_{c},\varphi(X)\rangle$ (14)

So inner loop in GBML can be interpreted as making proptotype vector with
given support set. $\square$

Table 8: Test accuracy % of 4-conv network on benchmark data sets. The values
in parentheses indicate the number of shots. The best accuracy among different
methods is bold-faced. To differentiate the notation, we have denoted
$\text{SHOT}_{3}$ as a model that uses 3 optimization steps in $f_{r}$ and 1
optimization step in $f_{t}$, and $\text{SHOT}_{6}$ as a model that uses 6
optimization steps in $f_{r}$ and 3 optimization steps in $f_{t}$.

meta-train miniImageNet Cars meta-test miniImageNet tieredImageNet Cars Cars
miniImageNet CUB MAML (1) 47.88 $\pm$ 0.55 51.84 $\pm$ 0.24 34.41 $\pm$ 0.47
47.78 $\pm$ 0.99 28.67 $\pm$ 1.17 30.95 $\pm$ 1.41 MAML + $\text{SHOT}_{3}$
(1) 47.97 $\pm$ 0.71 51.68 $\pm$ 0.68 34.03 $\pm$ 1.11 50.44 $\pm$ 0.62 30.19
$\pm$ 0.50 31.21 $\pm$ 0.64 MAML + $\text{SHOT}_{6}$ (1) 48.15 $\pm$ 0.31
51.67 $\pm$ 0.73 34.79 $\pm$ 0.70 49.89 $\pm$ 0.17 28.39 $\pm$ 0.37 30.83
$\pm$ 0.26 MAML (5) 64.81 $\pm$ 1.63 67.96 $\pm$ 1.22 46.57 $\pm$ 0.53 62.24
$\pm$ 2.01 37.23 $\pm$ 1.95 41.84 $\pm$ 1.25 MAML + $\text{SHOT}_{3}$ (5)
66.86 $\pm$ 0.58 69.48 $\pm$ 0.17 48.42 $\pm$ 0.72 69.08 $\pm$ 0.31 40.79
$\pm$ 0.93 43.46 $\pm$ 0.89 MAML + $\text{SHOT}_{6}$ (5) 66.27 $\pm$ 0.23
69.60 $\pm$ 0.31 45.83 $\pm$ 1.82 66.37 $\pm$ 1.93 39.35 $\pm$ 0.74 42.63
$\pm$ 0.24

##  Appendix E SHOT with more optimization step in the inner loop

In the main paper, we used only one step for the target model to improve
computation efficiency. However, it’s also important to test if SHOT works
with more optimization steps in the inner loop. As a reference model, we set
the number of optimization steps to 6, which is different from the main paper
where we only used 3 steps (same as the baseline). As shown in Table 9, SHOT
still performs better than the baseline even with more optimization steps in
the inner loop.

##  Appendix F SHOT can act as a regularlizer

Regularlization techniques are widely used in many GBML algorithms Rajeswaran
et al. [2019], Rusu et al. [2018a] and SHOT can also act as a regularizer. We
replaced the regularization term of LEO Rusu et al. [2018a], a popular GBML
algorithm that uses regularization in the inner loop, with SHOT. The results,
shown in Table 9, demonstrate that SHOT slightly enhances LEO’s performance.
Although the results are marginal, there is still room for improvement
considering LEO works with highly-fine-tuned hyperparameters. This suggests
the potential of SHOT as a plug-and-play regularization technique for GBML
algorithms, whether it is highly-fine-tuned or not. Also, the role of SHOT as
a regularizer explains why regularization techniques are effective in GBML. If
we minimize the total transport distance in the inner loop, parameters do not
move much even with a large learning rate, thereby the effect of the Hessian
reduces.

Table 9: Test accuracy % of LEO [Rusu et al., 2018a]. The values in parentheses are the number of shots. The better accuracy between the baseline and SHOT is bold-faced. Dataset | miniImageNet | tieredImageNet
---|---|---
LEO (1) | 60.24 $\pm$ 0.02 | 65.07 $\pm$ 0.14
LEO + $\text{SHOT}_{r}$ (1) |  60.39 $\pm$ 0.07 |  65.26 $\pm$ 0.10
LEO (5) | 75.27 $\pm$ 0.13 |  79.87 $\pm$ 0.06
LEO + $\text{SHOT}_{r}$ (5) |  75.32 $\pm$ 0.17 | 79.84 $\pm$ 0.14

##  Appendix G Algorithm

Algorithm 1 Training Algorithm for GBML Model

    Initialize parameter theta (the model)

    while not reach max number of epochs:
        # Sample query set Q, support set S
        # Start of Inner loop
        Initialize reference model and target model with theta
        theta_r = theta
        theta_t = theta

        # Optimize theta_r with more (e.g. 3) inner loop steps step_r
        # and a smaller learning rate alpha
        for i in range(step_r):
            theta_r = theta_r - alpha * gradient(loss(theta_r, S), theta_r)

        # Optimize theta_t with less (e.g. 1) inner loop steps step_t
        # and a larger learning rate step_r / step_t * alpha
        theta_t = theta_t - step_r / step_t * alpha *
            gradient(loss(theta_t, S), theta_t)
        # End of Inner loop

        # Start of Outer loop
        # Calculate outer loss
        outer_loss = loss(theta_r, Q)
        # Reduce the distance between theta_r and theta_t using KL-divergence
        kl_loss = KL_divergence(target_distribution(theta_t, Q)
            , reference_distribution(theta_r, Q).detach())
        outer_loss += lambda * kl_loss
        # Use meta-optimizer to optimize parameter theta
        theta = theta - beta * gradient(outer_loss, theta)
        # End of Outer loop
    # End of training loop

##  Appendix H Loss surface actually gets linear along the optimization
trajectory

Figure 4:  Comparison of loss surface between randomly initialized model and
$\text{SHOT}_{p}$ initialized model. Showing that effect of the Hessian is
suppressed is along the optimization trajectory.

Fig. 4 shows that our algorithm actually do what we expected. We compared
between randomly-initialized model and $\text{SHOT}_{p}$ initialized model.
Randomly-initialized model shows rough loss surface. The characteristic
doesn’t differ whether it is gradient’s direction or not. However,
$\text{SHOT}_{p}$ pretrained models shows smooth surface along the gradient
direction, while other direction remains rough.


# Reservoir Computing Using Complex Systems

N. Rasha Shanaz  Department of Physics, Bharathidasan University,
Tiruchirappalli 620024, India  K. Murali  Department of Physics, Anna
University, Chennai 600025, India  P. Muruganandam  Department of Physics,
Bharathidasan University, Tiruchirappalli 620024, India

###### Abstract

Reservoir Computing is an emerging machine learning framework which is a
versatile option for utilising physical systems for computation. In this
paper, we demonstrate how a single node reservoir, made of a simple electronic
circuit, can be employed for computation and explore the available options to
improve the computational capability of the physical reservoirs. We build a
reservoir computing system using a memristive chaotic oscillator as the
reservoir. We choose two of the available hyperparameters to find the optimal
working regime for the reservoir, resulting in two reservoir versions. We
compare the performance of both the reservoirs in a set of three non-temporal
tasks: approximating two non-chaotic polynomials and a chaotic trajectory of
the Lorenz time series. We also demonstrate how the dynamics of the physical
system plays a direct role in the reservoir’s hyperparameters and hence in the
reservoir’s prediction ability.

##  I Introduction

The emergence of chaos theory showed that the discrepancies in predicting and
forecasting highly nonlinear systems, like weather and economy, are due to the
chaos inherent in them. For years, this sensitivity to initial conditions
hindered long term predictions of complex systems. With the recent
developments in machine learning and artificial neural networks, forecasting
nonlinear chaotic data is becoming plausible and more accurate. These
computing paradigms show remarkable potential in predicting highly nonlinear
data, which was otherwise impossible.

One such paradigm, gaining momentum at present, is Reservoir Computing (RC)
[1]. The RC architecture consists of a reservoir, a randomly connected
artificial recurrent neural network (RNN), which projects the input into a
higher dimensional space and a readout layer, which analyses patterns in the
reservoir states. The reservoir is kept fixed, and only the readout is trained
[2]. This makes training, computationally inexpensive and the reservoir,
highly adaptable. The flexible architecture of RC enables a wide range of
physical systems to be used in place of RNNs as computing substrates [3].
Analog electronic circuits [4], discrete dynamical systems [5], optical
systems [6], and biological cells [7] are a few of the diverse examples. The
RC method has proven competence in forecasting spatiotemporally chaotic
dynamics, in the absence of defining models and only with observational data
for the system’s evolution [8].

The inherent dynamic nature of the RC paradigm renders itself an appropriate
candidate for intrinsically dynamic tasks. This paper aims to integrate the
model-free prediction ability of reservoir computing with the intrinsic
computational potential of complex dynamical systems. We computationally
implement a reservoir using a chaotic memristive circuit and study its
efficiency in approximating two nonlinear, non-chaotic polynomials and the
chaotic Lorenz time series. From the available hyperparameters, we built two
reservoirs and analysed how the corresponding dynamics play a role in the
efficiency of approximation. We find that the reservoir performs exceptionally
well in approximating non-chaotic data but forecasting the chaotic time series
requires drastic modification in the reservoir parameters.

###  I.1 Reservoir Computing architecture

The basic framework of a reservoir computing system consists of an input
layer, the reservoir, and the output or readout layer, as shown in Figure 1.

Figure 1: Schematic diagram of a general reservoir computing structure.

Conventionally, the reservoir is a randomly connected RNN. But any dynamical
system with a potential to dynamically respond to inputs can function as a
reservoir [3].

The reservoir converts the inputs given at the input layer into a
spatiotemporal pattern. The recurrent connections in the reservoir bring about
an intrinsic memory effect and nonlinearly map the input into a higher-
dimensional space. This transformation of input into higher-dimensional space
resembles what a complex dynamical system could do, hence introducing the
possibility of utilising complex systems themselves as reservoirs. The
reservoir states are the dynamical response of the reservoir system, which are
then mapped to the readout layer. The readout layer is trained using a
supervised, non-temporal mapping task like a linear regression [9]. The
reservoir offers a sufficiently rich cache of responses for an array of tasks.
This rules out the necessity to train the reservoir at every stage and for
each different task.

###  I.2 Optimal operating regime

For the reservoir system to efficiently perform computational tasks it should
have a high dimensionality for the mapping inputs with reservoir states, and
nonlinearity for extracting nonlinear dependencies of inputs. It should have
fading memory, also known as the echo state property, which means that the
influence of past inputs on the current reservoir states and outputs
asymptotically fades out [1, 10]. The mathematical definition says that a
system has echo state property if different initial states converge [11].
Along with these, it should have separation property, which maps different
inputs onto different reservoir states, and approximation property, which maps
infinitesimally different reservoir states onto identical targets.

Traditionally, this is realised through a recurrent network of numerous
reservoir nodes, but can also be achieved by time-multiplexing the output of a
single nonlinear node. This makes the state space infinite dimensional, but
the dynamics of the system remains finite dimensional [12].

From the complex system studies point of view, the reservoir can be regarded
as a complex dynamical system that operates optimally in a certain dynamical
regime. When a reservoir system functions in stable dynamical regions, small
variations in initial conditions will not reflect current dynamics, reducing
the complexity of the model and increasing the bias. This might result in
underfitting of data. Whereas, when it functions in chaotic regimes,
sensitivity increases and amplifies infinitesimal variations, compromising the
approximation of closer values and hence increasing variance. This could cause
overfitting and the reservoir will not be versatile for different data sets.
If the reservoir is operated in a dynamical region, where the dynamics are
neither too contracting nor too diverging, it will be able to pick up
variations in data and also bring about a sweeping approximation of
hypotheses. Such an optimal working regime can be obtained at the edge of
chaos, a region where both order and chaos coexist. [13, 14].

In the wide variety of available nonlinear dynamical systems, the motive is to
choose a simple system that is capable of exhibiting rich nonlinear dynamics
and has options for more hyperparameters allowing more precision. A system
with chaotic dynamics is an ideal choice for a reservoir because it easily
fulfils the above mentioned criteria.

##  II Methods

###  II.1 The reservoir

Based on the above mentioned properties of the reservoir, we chose a chaotic
oscillator circuit with a memristive element.

Figure 2: The second order memristive oscillator used as the reservoir

It is a second order non-autonomous memristive oscillator (Figure 2), built by
modifying the memristive Chua’s circuit [15]. We took into consideration, its
simple circuit architecture, experimental reproducibility, rich dynamics,and
the inherent fading memory property of the memristor [16].

Figure 3: Bifurcation diagrams of $v_{1}$ with (a) increasing resistance $R$
and (b) increasing forcing amplitude $A$.

The governing equations of the circuit are:

$\displaystyle\dfrac{dv_{1}}{dt}$
$\displaystyle=\dfrac{v_{s}-v_{1}}{RC_{1}}+\dfrac{(1-gv_{0}^{2})v_{1}}{R_{3}C_{1}},$
(1a) $\displaystyle\dfrac{dv_{0}}{dt}$
$\displaystyle=-\dfrac{v_{1}}{R_{1}C_{0}}-\dfrac{v_{0}}{R_{2}C_{0}},$ (1b)

where $v_{s}=A\sin(2\pi\nu t)$ is the forcing signal with amplitude $A$ and
angular frequency $\omega=2\pi\nu$. The memristor ($W)$ is modeled based on
[17] as

$\displaystyle i=W(v_{0})v=-\dfrac{1}{R_{3}}\left(1-gv_{0}^{2}\right)v,$ (2)

where $W(v_{0})$, $i$ and $v$ are the memductance, current and voltage of the
memristor, respectively, and $g$ is the internal gain factor. For
computational convenience, we use the following set of normalized version of
the governing equations.

$\displaystyle\dot{x}$ $\displaystyle=-\alpha y-x,$ (3a)
$\displaystyle\dot{y}$
$\displaystyle=-\beta\left(A^{\prime}\sin\omega^{\prime}t-y\right)+\gamma\left(1-x^{2}\right)y,$
(3b)

where $v_{0}=x/\sqrt{g}$, $v_{1}=y/\sqrt{g}$, $A={A^{\prime}}/{\sqrt{g}}$,
$t=R_{2}C_{0}t^{\prime}$, $\alpha={R_{2}}/{R_{1}}$,
$\beta={R_{2}C_{0}}/{RC_{1}}$, $\gamma={R_{2}C_{0}}/{(R_{3}C_{1})}$, and
$\omega={\omega^{\prime}}/{(R_{2}C_{0})}$. Here, $x$ corresponds to the
dimensionless voltage across the capacitor $C$ and $y$ corresponds to the
dimensionless voltage through the memristor $W$. The term $A\sin\omega t$ is
the external forcing signal with amplitude $A$ and angular frequency $\omega$.
The values of the memristor parameters were set as
$R_{1}=8.0~{}\mbox{k}\Omega$, $R_{2}=4.0~{}\mbox{k}\Omega$,
$R_{3}=1.4~{}\mbox{k}\Omega$, $C_{0}=4.7~{}nF$, $C_{1}=6.8~{}nF$ and $g=0.1$
[15].

The dynamics of the circuit is governed by forcing amplitude ($A$), resistance
($R$), and forcing frequency ($\nu$). The availability of more parameters is
advantageous for fine tuning the reservoir computing abilities of the circuit.
Keeping in mind the experimental hardships, we chose $R$ and $A$ to be the
hyperparameters for the RC. Figure 3(a) shows the bifurcation diagrams of
$v_{1}$ with $R$ increasing from $1.9~{}\mbox{k}\Omega$ to
$2.8~{}\mbox{k}\Omega$, with $A$ fixed at $2~{}V$, and figure 3(b) shows the
bifurcation of $v_{1}$ with $A$ increasing from $1.5~{}V$ to $3.5~{}V$,
keeping $R$ fixed at $2.6~{}\mbox{k}\Omega$.

A wide range of striking dynamical behaviour can be observed, including
periodic limit cycles, period doubling bifurcations and chaos. A few selected
dynamics are presented by the numerically simulated trajectories of (1) on the
$v_{0}-v_{1}$ plane, for arbitrary values of $R$, as shown in figure 4.

Figure 4: The numerical phase portraits for different values of $R$ on the
$v_{0}-v_{1}$ plane. (a) Period-$T$ limit cycle at $R=1.9\mbox{k}\Omega$. (b)
Chaotic attractor at $R=2.1\mbox{k}\Omega$. (c) Period-$3T$ limit cycle at
$R=2.3\mbox{k}\Omega$. (d) Chaotic spiral attractor at $R=2.7\mbox{k}\Omega$.

The dynamics corresponding to $A$ have similar features. Since such
environment paves way for the reservoir to be operated in the edge of chaos,
the existence of chaotic and periodic windows is essential for this system to
be employed as a reservoir.

To analyse the effect of the two hyperparameters on the computational ability
of RC, we built two reservoirs: Resistance controlled ($R$-Reservoir) and
Amplitude controlled ($A$-Reservoir) versions. For nonchaotic, nonlinear
approximating task, we chose the working regime of the $R$-reservoir to be 2.2
$\mbox{k}\Omega-2.8~{}\mbox{k}\Omega$ and that of $A$-reservoir to be
$2.1V-2.5V$, which are regions with distinct yet coexisting chaotic and stable
dynamics. And for the case of approximating chaotic data, we operated the
$R$-reservoir between $2.10~{}\mbox{k}\Omega-2.19~{}\mbox{k}\Omega$, and the
$A$-reservoir between $2.75V-3.25V$, which are regions with predominantly
chaotic dynamics. The reasoning behind this choice will be explained in
section III

###  II.2 The task

We evaluate the performance of the reservoir on a non-temporal task of
approximating three different data sets of varying complexity. We estimate the
outcome for a new data set based on the training on a known data set. The
three data sets used for the task are:

  * •

$5^{th}$-degree polynomial:

$\displaystyle f(x)=x^{5}-5x^{4}+5x^{3}+5x^{2}-6x-1$ (4)

with $x$ ranging from $-1.25$ to $3.25$

  * •

$9^{th}$-degree polynomial:

$\displaystyle f(x)=$ $\displaystyle\,x^{9}+3x^{8}-4.5x^{7}-21x^{6}+x^{5}$
$\displaystyle\,+44x^{4}+13x^{3}-25x^{2}-11x,$ (5)

with $x$ ranging from $-1.791$ to $1.834$

  * •

The Lorenz time series [18]:

$\displaystyle\dot{x}$ $\displaystyle=\sigma(y-x),$ (6a)
$\displaystyle\dot{y}$ $\displaystyle=\rho x-y-xz,$ (6b)
$\displaystyle\dot{z}$ $\displaystyle=xy-\beta z,$ (6c)

where $\sigma=10$, $\rho=28$ and $\beta=2.667$, with the initial condition set
at $x_{0}=0.5$, $y_{0}=1$, and $z_{0}=2$. We approximate 10,000 data points of
the trajectory $x(t)$, leaving sufficient points as transients.

###  II.3 Mapping

The prime concern is to define computation in terms of dynamics of the system.
We employ a general framework based on the reservoir computing prediction
methods proposed by Jaeger and Haas [2]. For an input vector $u(t)$, the
reservoir state vector is modelled as $r(t)=F[W_{in}~{}u(t)+M\,r(t-1)]$, where
$t$ denotes discrete time, $W_{in}$ is the input weight matrix, $M$ is the
reservoir weight matrix and $F$ is the activation function, which can be any
differentiable function and its form may not necessarily be known [19]. Here
the discrete-time parameter $t$ does not represent time in a real sense;
instead, it represents the parameter on which the input vector $u(t)$ depends.
The above equation represents a non-autonomous dynamical system driven by an
external input $u(t)$ [3]. The output vector is given by $y(t)=W_{out}r(t)$,
where, $W_{out}$ is the weight matrix of the readout, which is trained using
ridge regression. The reservoir is considered as a black box model of the
system employed [2].

We replace the standard neuron models for $F$,with the circuit’s equations
(3). We apply the input through the parameter $\beta$ for the $R$-Reservoir
and the external forcing signal $A\sin\omega t$ for the $A$-reservoir. The
reservoir input parameter $X$ is normalised and mapped to the input range $x$
of the polynomial. Here $X$ can be $R$ or $A$ in the respective reservoir
versions and is taken in a linear range $X_{min}$ to $X_{max}$, which are the
boundary values of the chosen working regime, as specified in section II.1.
For instance, for the $R$-reservoir in the task of approximating the
$5^{th}$-degree polynomial (4), the input is normalised by mapping
$[R_{min},R_{max}]$ to $u(t)$ as
$[2.1\,\mbox{k}\Omega,2.19\,\mbox{k}\Omega]\rightarrow[0,1]$. This normalised
input $u(t)$ is then mapped to the independent variable, $x$, of the
polynomial as $[0,1]\rightarrow[-1.25,3.25]$. The interval is divided into
$10000$ equally spaced points.

For each input $u$, the reservoir state vector is perturbed for $N=5$ periods
and for each cycle of the perturbing signal, $k=10$ points are chosen, to
generate $kN$ discrete outputs. This gives $kN=50$ outputs for every parameter
value.

Figure 5: The $kN=50$ points (blue circles) chosen from the reservoir system
output for $N=5$ periods of forcing signal. This output is taken for
$R=2.2~{}\mbox{k}\Omega.$ Figure 6: The approximated polynomials by the
$R$-reservoir (above, in red) and $A$-reservoir (below, in blue). The target
data is shown in cyan. Approximation of (a) the $5^{th}$-degree polynomial (4)
(MSE=$0.0766$), (b) the $9^{th}$-degree polynomial (5) (MSE=$0.2499$) and (c)
the Lorenz time series (6 (MSE=$11.1691$), by the $R$-reservoir. Approximation
of (d) the $5^{th}$-degree polynomial (MSE=$0.000102$), (e) the
$9^{th}$-degree polynomial (MSE=$0.001170$) and d(f) the Lorenz time series
(MSE=$9.0455$) by the $A$-reservoir.

The system is reset to start in the same initial condition between each input.
A visual representation of how the $kN$ points are chosen is shown in Figure
5. Since these points are outputs of a chaotic system, there will be an
inherent complexity in them. This brings about the mapping of input to a
higher dimensional space, by generating numerous nonlinear transformations as
a function of the input parameter $X$, without the requirement for a system
dynamics with higher dimensionality.

We use the Python based Scikit-learn package for ridge regression to train the
reservoir for these tasks [20]. We assess the performance of the reservoir
using Mean Squared Error (MSE) from cross-validation. The details are outlined
in A We optimise numerically using the train and test scheme to test the
goodness of the reservoir. We split the data set in a training:testing ratio
of $67:33$. In order for the algorithm to learn the overall pattern of the
data, the training points are chosen randomly from the entire data range. And
the remaining data points are used for the approximation task. We also study
the effect of the working parameter range of the reservoir on its bias and
variance and thereby analysing the best working regime for approximation
tasks.

##  III Results and discussion

The three data sets were approximated by the reservoir and the approximated
curves of the $5^{th}$\- and $9^{th}$-degree polynomials, and the Lorenz time
series are shown in Figure 6. All points shown on the plot are the
approximated values and the training data are not shown here. For comparison,
the target data is presented. For all three cases, the reservoir system
parameters were kept fixed. The goodness of fit is assessed through the loss
function, Mean Squared Error (MSE).

Figure 7: The various approximations of the Lorenz time series by the
reservoirs in different dynamical regions, demonstrating the effect of the
system dynamics on computational capability of the reservoir. The plots show
approximation by the $R$–reservoir (above, in red) operated in (a) periodic
regime ($R=2.3\mbox{k}\Omega-2.4\mbox{k}\Omega$ ) (b), chaotic regime
$R=2.10\mbox{k}\Omega-2.19\mbox{k}\Omega$ and (c) the full bifurcation
parameter range ($R=1.9\mbox{k}\Omega-2.8\mbox{k}\Omega$), and the
approximation by the $A$–reservoir (below, in blue) in (d) periodic regime
($A=1.5V-1.8V$), (e) chaotic regime ($A=2.5V-3.4V$) and (f) the full
bifurcation parameter range ($A=1.5V-3.5V$).

The variances of the reservoirs were regulated to minimise the loss function
and to accommodate the increasing complexity in the three data sets and
thereby, minimising error.

The performance of the reservoir, in terms of approximation, is good in the
case of the $5^{th}$\- and $9^{th}$-degree polynomials. Comparing the mean
square errors for the two reservoirs, the $A$-reservoir has a very low value
of the order $10^{-4}$ for $5^{th}$-degree and $10^{-3}$ for the
$9^{th}$-degree polynomial. This is two orders of magnitude less than the MSEs
of the $R$-reservoir, which are $10^{-2}$ and $10^{-1}$ for the $5^{th}$\- and
$9^{th}$-degree polynomials. The absence of noise in the approximation and the
low values of MSE show that the $A$-reservoir is better performing. This can
be attributed to the fact that the dynamical regime, in which the A-reservoir
was operated in ($2.0V-2.5V$) has a good balance of order and chaos, and thus
bringing the optimal balance between bias and variance.

Now, in the case of Lorenz time series, we plot 10,000 data points of the time
series for $10$ units of time. In the approximation tasks, the MSEs are too
high but the overall shape of the approximated polynomial resembles the time
series, for both the reservoirs. Improving the quality of approximation, in
this case, will require increasing the variance to very high values in order
to attune to the chaos. Such an increment will result in the reservoir
capturing the minute variations along with the underlying pattern. This will
lead to overfitting of other data sets, hence compromising the long term
prediction ability of the reservoir. This led us to study the effect of
different dynamical regions on the computational ability and thereby eliciting
the optimal working regime for predicting chaotic data.

Taking the same time series, we trained both the A-controlled and
$R$-controlled reservoirs in three different dynamical regimes of the
memristive circuit. We chose the regions with only periodic dynamics, only
chaotic dynamics, and the entire bifurcation range. Figure 7 shows, the
various approximations by the $R$-reservoir and the $A$-reservoir and the
corresponding parameter ranges in which they were operated.

The reservoir system dynamics for a particular range of hyperparameters ($R$
and $A$), as seen in the bifurcation diagram, is reflected in the
corresponding approximation of the Lorenz equation. For instance, for the
Figure 7(c) the reservoir parameter $R$ is set in the range
$1.9\mbox{k}\Omega-2.8\mbox{k}\Omega$. The bifurcation diagram of $R$ given in
Figure 3(a) shows alternating periodic and chaotic regions. In the figure, we
find that the approximated curve is flatly approximated in ranges
corresponding to periodic dynamics and is highly sensitive in ranges
corresponding to chaotic dynamics, mirroring the system dynamics. Similar
interpretations can be drawn by comparing the bifurcation plot region for a
particular parameter range with the corresponding plot of approximation done
by the reservoir operating in that range.

This demonstrates how system dynamics will alter bias and variance of a
reservoir built of a complex system. In the stable regime, a system will
behave predictably and has a good fading memory. Small variations in initial
conditions or noise will not be reflected in the current dynamics. Hence along
with rejecting noise, a reservoir operating in this region will fail to
capture the underlying pattern in the data. On contrary, a system in the
chaotic regime will have high sensitivity to infinitesimal changes and a
diverging behaviour. This will compromise the reservoir’s ability to
approximate close-lying values and hence rendering it incapable of predicting
other data sets. And hence a system with highly chaotic dynamics will be
better able to predict a chaotic data set.

##  IV Conclusion

We demonstrated how a single node reservoir, made of a simple memristive
chaotic oscillator, can be employed to approximate higher degree polynomials
and how its dynamics play a role in the bias and variance of the reservoir’s
computational abilities. Of the two reservoirs built, the amplitude controlled
version performs better due to the optimal working region offered by the
parameter. The approximated data matches the target data set with a very small
mean squared error. By operating the reservoir in different dynamical regimes,
we re-established how chaos and order affect the variance and bias of the
reservoir.

This work is limited to non-temporal tasks but holds open the future
possibility for temporal time series forecasting. The presence of memristor in
the circuit makes it a promising candidate for human neural analysis and
studying complex adaptive systems.It is worthy to note that reservoirs made of
memristor networks satisfy fading memory and separation properties [21].

There are many real-world scenarios where order and chaos play an important
role, ranging from epileptic seizures to binary star systems. The inherent
nonlinearity in these systems makes long term prediction difficult. The
current scenario in machine learning and artificial neural networks poses a
promising future. Incorporating the dynamics of nonlinear complex systems will
take forward the current developments by many steps.

##  Appendix A Regression and cross validation

###  A.1 Ridge regression

In machine learning, regression analysis is performed to estimate the
relationship between the independent variable and dependant outcome
variable(s). It establishes a relationship between the input and output
variables of the reservoir and hence helps in prediction or forecasting tasks.
The most common method is linear regression where we find a line that most
closely fits the data. For a linear model of $i$ number of features $X_{i}$
given as $Y_{i}=\underset{i}{\sum}w_{i}X_{i}+b$, the regression aims to
optimise the slope $w_{i}$ and the intercept $b$ by minimising the cost
function given as

$\displaystyle\underset{w}{\text{min\
}}\left(\underset{i}{\sum}(w_{i}X_{i}-Y_{i})^{2}\right).$ (7)

Ridge regression is a type of linear regression which uses regularisation and
seeks to minimise the sum of squared error between the model and the training
data. The ridge coefficients minimise the cost function by imposing a penalty
on the coefficient sizes equal to the square of their magnitudes as

$\displaystyle\underset{w}{\text{min\
}}\left(\underset{i}{\sum}\left(w_{i}X_{i}-Y_{i}\right)^{2}+\alpha~{}\underset{i}{\sum}w_{i}^{2}\right),$
(8)

where $\alpha\geq 0$ is the complexity parameter.

###  A.2 Cross-validation

It is a technique for assessing the validity of the model learned by the
machine learning algorithm. It is done by partitioning the data into
complementary subsets, analysing one _training_ subset and validating the
analysis on the other _validation_ subset. Performing this multiple times will
help reduce the variance. The cross-validation estimator will select the best
hyperparameter value.

## Acknowledgement

The work of PM forms parts of research projects by Council of Scientific and
Industrial Research (CSIR), India under Grant No. 03(1422)/18/EMR-II, and
Science and Engineering Research Board (SERB), India under Grant No.
CRG/2019/004059.

## References

  * [1] H Jaeger,  GMD Report 148, 1 (2001).
  * [2] H Jaeger and H Haas,  Science 304, 78 (2004).
  * [3] G Tanaka, T Yamane, J B Héroux, R Nakane, N Kanazawa, S Takeda, H Numata, D Nakano, and A Hirose,  Neural Networks 115, 100 (2019).
  * [4] J H Jensen and G Tufte,  Artificial Life Conference Proceedings 29, 222 (2017).
  * [5] D Snyder, A Goudarzi, and C Teuscher,  Phys. Rev. E 87, 042808, 2013.
  * [6] K Vandoorne, W Dierckx, B Schrauwen, D Verstraeten, R Baets, P Bienstman, and J V Campenhout,  Optics express 16, 15 (2008).
  * [7] B Jones, D Stekel, J Rowe, and C Fernando,  Proc. IEEE Symposium on Artificial Life,187 (2007).
  * [8] J Pathak, B Hunt, M Girvan, Z Lu, and E Ott,  Phys. Rev. Lett. 120 (2018).
  * [9] M Lukoševičius and H Jaeger,  Comput. Sci. Rev. 3, 127 (2009).
  * [10] W Maass, T Natschläger, and H Markram,  Neural Comput. 14, 2531 (2002).
  * [11] J Boedecker, O Obst, JT Lizier, N M Mayer, M Asada,  Theory in Biosciences 131, (2012).
  * [12] L Appeltant, M C Soriano, G Van Der Sande, J Danckaert, S Massar, J Dambre, B Schrauwen, C R Mirasso, and I Fischer,  Nat. Commun. 2 (2011).
  * [13] R Legenstein and W Maass,  Neural Networks 20, 323 (2007).
  * [14] M Inubushi and K Yoshimura,  Sci. Rep. 7, 1 (2017).
  * [15] Q Xu, Q Zhang, B Bao, and Y Hu,  IEEE Access 5, 21039 (2017).
  * [16] A Ascoli, R Tetzlaff, L O Chua, J P Strachan, and R S Williams,  IEEE Trans. Circuits Syst. I, Reg. Papers 63, 389 (2016).
  * [17] Q Xu, Y Lin, B Bao, and M Chen,  Chaos Soliton Fract. 83, 186 (2016).
  * [18] E H Lorenz,  J. Atmos. Sci. 20, 131, (1963).
  * [19] Z Lu, B R Hunt, and E Ott,  Chaos 28, 061104 (2018).
  * [20] F Pedregosa et al.,  _J. Mach. Learn. Res._ 12, 2825 (2011).
  * [21] F Sheldon, A Kolchinsky and F Caravelli _Phys. Rev. E_ 106, 045310 (2022).