# Solar cycle 25 prediction using length-to-amplitude relations

Vladimir G. Ivanov1,2 1Central Astronomical Observatory at Pulkovo, Saint-
Petersburg, Russia
2email<EMAIL_ADDRESS>
###### Abstract

We propose a simple method for prediction of the 11-year solar cycle maximum
that is based on two relations. One of them is well known Waldmeier’s rule
that binds the amplitude of a cycle and the length of its ascending phase. The
second rule relates the length of a given cycle from minimum to minimum and
the amplitude of the next one. Using corresponding linear regressions we
obtain for the amplitide of cycle 25 in the scale of 13-month smoothed monthly
total revised sunspot number ${\rm SN}_{\rm max}(25)=181\pm 46$ and for the
moment of the maximum $T_{\rm max}(25)=2024.2\pm 1.0$. Therefore, according to
the prediction, cycle 25 will be higher than the previous one (${\rm SN}_{\rm
max}(24)=116$) with probability 0.92.

###### keywords:

solar cycle prediction; cycle 25; Waldmeier rule

[Du et al. (2022)] du22

##  1 Introduction

In 2019 experts of Solar Cycle 25 Prediction Panel predicted that the maximum
of the solar cycle 25 (SC25) will lie in the sunspot number range between 95
and 130, that is similar to SC24 with the maximal SN=116
(https://www.weather.gov/news/190504-sun-activity-in-solar-cycle). The opinion
that SC25 will be equal or lower than SC24 is shared by many authors (du22;
burud21; nandy21; wu21; courtillot21 etc).

If the odd SC25 is lower than the even SC24, it will violate the Gnevyshev-Ohl
correlation rule [Gnevyshev and Ohl (1948)]. Many regard this as a sign that
the Sun is entering to a Dalton-like global minimum.

However, there is an alternative point of view that the forecast of the
experts has been underestimated and the forecoming cycle will be higher
(macintosh20; koutchmy21; prasad22; lu22 etc). None of the two positions have
prevailed so far.

In the following we propose a simple method of solar maximums predictions
based on linear relations between amplitudes of cycles and lengths of its
phases. This method provides one more argument that SC25 will be higher that
SC24 and comparable with SC23.

##  2 Data and notation

Hereafter we will use for analysys and prediction the 13-month smoothed
monthly averages of the recalibrated sunspot number SN for 1749–2021 [Clette
et al. (2014)] (https://wwwbis.sidc.be/silso/DATA/SN_ms_tot_V2.0.txt).

We will introduce the followig notation: ${\rm SN}_{\rm max}(i)$ is the
sunspot index in the maximum of the $i$th cycle (i.e. its amplitude), $T_{\rm
min}(i)$ and $T_{\rm Max}(i)$ are moments of its minimum and maximum, $T_{\rm
mm}(i)=T_{\rm min}(i+1)-T_{\rm min}(i)$ is its length (from minimum to
minimum), $T_{\rm mM}(i)=T_{\rm Max}(i)-T_{\rm min}(i)$ is the length of the
ascending branch of the cycle (from minimum to maximum).

##  3 The Waldmeiers rule

The well-known empirical Waldmeier rule (WR) states that the length of the
ascending branch of the cycle anticorrelates with its amplitude [Waldmeier
(1935)] (Fig. 1). One can find the linear regression equation

$\displaystyle{\rm SN}_{\rm max}(i)$ $\displaystyle=$ $\displaystyle
a_{1}+b_{1}\,T_{\rm mM}(i)\,,$ $\displaystyle a_{1}$ $\displaystyle=$
$\displaystyle 340.2\pm 31.3\,,$ $\displaystyle b_{1}$ $\displaystyle=$
$\displaystyle(-36.76\pm 6.90)\;{\rm yr}^{-1}\,,$ (1)

and the correlation coefficient between the regression parameters is

$r(a_{1},b_{1})=-0.976\,.$

Figure 1: The relationship between the length of the ascending phase $T_{\rm
mM}(i)$ and the amplitude of the cycle ${\rm SN}_{\rm max}(i)$ (the Waldmeier
rule). The correlation coefficient is $r=-0.75$. The solid curve is the linear
regression (1). The dashed curves mark the standard error ranges for
predictions.

##  4 The length-to-next-amplitude rule (LNAR)

Another empirical rule, similar to WR but less known (and still nameless, to
the best of our knowledge), states that the length of a given cycle
anticorrelates with the amplitude of the next one [Hathaway et al. (1994),
Solanki et al. (2002), Hathaway (2015), Ivanov (2021)] (Fig. 2). The
corresponding regression is

$\displaystyle{\rm SN}_{\rm max}(i+1)$ $\displaystyle=$ $\displaystyle
a_{2}+b_{2}\,T_{\rm mm}(i)\,,$ $\displaystyle a_{2}$ $\displaystyle=$
$\displaystyle 544.5\pm 86.6\,,$ $\displaystyle b_{2}$ $\displaystyle=$
$\displaystyle(-33.01\pm 7.81)\;{\rm yr}^{-1}\,,$ $\displaystyle
r(a_{2},b_{2})$ $\displaystyle=$ $\displaystyle-0.994\,.$ (2)

Figure 2: The relationship between the length of the cycle $T_{\rm mm}(i)$ and
the amplitude of the next one ${\rm SN}_{\rm max}(i+1)$ (the LNA rule). The
correlation coefficient is $r=-0.68$. The solid curve is the linear regression
(2). The red circle is the prediction for SC25. The dashed curves mark the
standard error ranges for predictions.

##  5 The length-to-ascending-length rule (LALR)

A direct consequence of the WR and LNAR is the third rule, that binds the
length of the current cycle with the length of the ascending phase of the next
one (Fig. 3). These parameters correlate, and the regression is

$\displaystyle{\rm T}_{\rm mM}(i+1)$ $\displaystyle=$ $\displaystyle
a_{3}+b_{3}\,T_{\rm mm}(i)\,,$ $\displaystyle a_{3}$ $\displaystyle=$
$\displaystyle(-1.98\pm 1.81)\;{\rm yr}\,,$ $\displaystyle b_{3}$
$\displaystyle=$ $\displaystyle 0.570\pm 0.163\,,$ $\displaystyle
r(a_{3},b_{3})$ $\displaystyle=$ $\displaystyle-0.994\,.$ (3)

Figure 3: The relationship between the length of the cycle $T_{\rm mm}(i)$ and
the length of the next ascending phase $T_{\rm mM}(i+1)$ (the LAL rule). The
correlation coefficient is $r=+0.61$. The red circle is the prediction for
SC25. The solid curve is the linear regression (3). The dashed curves mark the
standard error ranges for predictions.

##  6 Prediction of SC25

The LNA and LAL rules stated above allow to estimate the time and magnitude of
the maximum of the next cycle provided we know the length of the current one.
Since the minimum between SC24 and SC25 in the 13-month smoothed index occurs
in December 2019, the length of the SC24 is 11.0 yr, and it follows from (1)
and (3) that

${\rm SN}_{\rm max}(25)=181\pm 46\,,$ ${\rm T}_{\rm mM}(25)=(4.23\pm
0.96)\;{\rm yr}\,,$

and

$T_{\rm max}(25)=(2024.24\pm 0.96)\;{\rm yr}$

(see Figs. 2 and 3).

The errors of the predictions are large enough. Nevertheless, according to the
prediction, with probability 0.92 SC25 will be higher than SC24 (${\rm
SN}_{\rm max}(24)=116$).

Probably, it will be similar to SC23 (${\rm SN}_{\rm max}(23)=176$) or SC17
(${\rm SN}_{\rm max}(17)=189$) (Fig. 4). It means that the odd SC25 will be
higher than the even SC24, and the Gnevyshev-Ohl rule will keep valid.

Figure 4:  The observed SN index for SC24 and SC25 (black line), SC17 (short
dashes) and SC23 (long dashes) shifted to match the minimums to that of SC25,
our prediction of SC25 maximum (red point with errorbars), the prediction of
SC25 (green line) and its maximum (green point with errorbars) by the Solar
Cycle Prediction Panel (https://services.swpc.noaa.gov/json/solar-
cycle/predicted-solar-cycle.json).

Earlier [Ivanov (2021)], using the same method, we have obtained for SC25
weaker estimates: the amplitude $136\pm 36$ and the moment of maximum
$2025.7\pm 0.7$ (for the Gaussian smoothing of SN with $\sigma=8$ months).
This underestimate can be explained by the fact that in January 2021, when the
mentioned article was being prepared, it was difficult to determine the moment
of SC25 minimum in the smoothed index accurately, and we assumed it to occur
too late, in October 2020, when the last local minimum of the index before its
fast growth took place.

##  7 Control of stability

To control the stability of the prediction method we will do the following
procedure. Let’s take the subseries of parameters from 1st to $(i-1)$th cycle,
construct the regressions for the LNAR and LALR using the truncated series and
obtain the predictions for the moment and amplitude of the $i$th cycle. The
results of the predictions compared with the observed values are plotted in
Figs. 5 and 6. All predicted values are in $1.5\sigma$ bands relative to the
observed values. The number of predictions in the one sigma band is 8 out of
13 (62%) for amplitudes and 11 out of 13 (85%) for lengths, which does not
contradict to the probability to be in the $1\sigma$ band for the normal
distribution (68%).

##  8 Conclusions

The proposed method of prediction, despite its simplicity and rather large
errors, provides an accuracy of about 15–30% for the moment of maximum and
25–50% in its magnitude, which is enough to distinguish between cycles of low,
medium and high power. In particular, the method predicts SC25 of medium
magnitude ${\rm SN}_{\rm max}(25)=181\pm 46$, which will be higher than the
previous one (${\rm SN}_{\rm max}(24)=116$) with probability 0.92, and with
the same probability the Gnevyshev-Ohl correlation rule will be valid for
SC25.

Figure 5: The control predictions for ${\rm SN}_{\rm max}$ on truncated
series. The red circles with errorbars are predictions of the $i$th SC. The
blue circles are observed values. Figure 6: The control predictions for
$T_{\rm mM}$ on truncated series. The red circles with errorbars are
predictions of the $i$th SC. The blue circles are observed values.

## References

  * Burud et al. (2021) Burud, D.S., Jain, R,, Awasthi, A.K., Chaudhari, S., Tripathy, S.C., Gopalswamy, N., Chamadia, P., Kaushik, S,C., Vhatkar, R.: 2021, Research in Astronomy and Astrophysics 21, 215.
  * Clette et al. (2014) Clette, F., Svalgaard, L., Vaquero, J.M., and Cliver, E.W.: 2014, Space Sci. Rev. 186, 35.
  * Courtillot et al. (2021) Courtillot, V., Lopes, F., Le Mouël, J.L.: 2021, Sol. Phys. 296, 21.
  * Du et al. (2022) Du, Z.L.: 2022 Ap&SS 367, 20.
  * Gnevyshev and Ohl (1948) Gnevyshev, M.N. and Ohl, A.I.: 1948, Astron. Zh. 25, 18.
  * Hathaway (2015) Hathaway, D.H.: 2015 Living Rev. Sol. Phys. 12, 4.
  * Hathaway et al. (1994) Hathaway, D.H., Wilson, R.M., and Reichmann, E.J.: 1994, Sol. Phys. 151, 177\.
  * Ivanov (2021) Ivanov, V.: 2021, Geomagnetism and Aeronomy 61, 1029.
  * Koutchmy et al. (2021) Koutchmy, S., Tavabi, E., Noëns, J.-C., Wurmser, O., Filippov, B.: 2021, In: Siebert A. et al. (eds.) Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics, 238.
  * Lu et al. (2022) Lu, J.Y., Xiong, Y.T., Zhao, K., Wang, V, Li, J.Y., Peng, G.S., and Sun, M.: 2022, ApJ 924, 59.
  * McIntosh et al. (2020) McIntosh, S.W., Chapman, S.C., Leamon, R.J., Egeland, R., Watkins, N.W.: 2020, Sol. Phys., 295, 163.
  * Nandy (2021) Nandy, D.: 2021, Sol. Phys. 296, 54.
  * Solanki et al. (2002) Solanki, S.K., Krivova, N.A., Schüssler, M., and Fligge, M.: 2002 A&A 396, 1029.
  * Prasad et al. (2022) Prasad, A., Roy, S., Sarkar, A., Chandra Panja, S., Narayan Patra, S.: 2022,  Adv. Space Res. 69, 798.
  * Waldmeier (1935) Waldmeier, M.: 1935, Astron. Mitt. Eidgenössischen Sternwarte Zürich 14, 105.
  * Wu and Qin (2021) Wu, S.-S., Qin, G.: 2021, E-print arXiv:2102.06001.


* [27] G. Ballesteros, J. Rey, M. Taoso, and A. Urbano, “Stochastic inflationary dynamics beyond slow-roll and consequences for primordial black hole formation,” JCAP 08 (2020) 043, arXiv:2006.14597 [astro-ph.CO].
  * [28] C. T. Byrnes, P. S. Cole, and S. P. Patil, “Steepest growth of the power spectrum and primordial black holes,” JCAP 06 (2019) 028, arXiv:1811.11158 [astro-ph.CO].
  * [29] S. Garcia-Saenz and S. Renaux-Petel, “Flattened non-Gaussianities from the effective field theory of inflation with imaginary speed of sound,” JCAP 11 (2018) 005, arXiv:1805.12563 [hep-th].
  * [30] T. Bjorkmo, R. Z. Ferreira, and M. D. Marsh, “Mild Non-Gaussianities under Perturbative Control from Rapid-Turn Inflation Models,” JCAP 12 (2019) 036, arXiv:1908.11316 [hep-th].
  * [31] J. Fumagalli, S. Renaux-Petel, J. W. Ronayne, and L. T. Witkowski, “Turning in the landscape: a new mechanism for generating Primordial Black Holes,” arXiv:2004.08369 [hep-th].
  * [32] R. Z. Ferreira, “Non-Gaussianities in models of inflation with large and negative entropic masses,” JCAP 08 (2020) 034, arXiv:2003.13410 [astro-ph.CO].
  * [33] L. Senatore, “Tilted ghost inflation,” Phys. Rev. D71 (2005) 043512, arXiv:astro-ph/0406187 [astro-ph].
  * [34] G. Panagopoulos and E. Silverstein, “Primordial Black Holes from non-Gaussian tails,” arXiv:1906.02827 [hep-th].
  * [35] J. M. Ezquiaga, J. García-Bellido, and V. Vennin, “The exponential tail of inflationary fluctuations: consequences for primordial black holes,” JCAP 03 (2020) 029, arXiv:1912.05399 [astro-ph.CO].
  * [36] D. G. Figueroa, S. Raatikainen, S. Rasanen, and E. Tomberg, “Non-Gaussian tail of the curvature perturbation in stochastic ultra-slow-roll inflation: implications for primordial black hole production,” arXiv:2012.06551 [astro-ph.CO].
  * [37] C. Pattison, V. Vennin, D. Wands, and H. Assadullahi, “Ultra-slow-roll inflation with quantum diffusion,” JCAP 04 (2021) 080, arXiv:2101.05741 [astro-ph.CO].
  * [38] M. Celoria, P. Creminelli, G. Tambalo, and V. Yingcharoenrat, “Beyond Perturbation Theory in Inflation,” arXiv:2103.09244 [hep-th].
  * [39] M. Biagetti, V. De Luca, G. Franciolini, A. Kehagias, and A. Riotto, “The Formation Probability of Primordial Black Holes,” arXiv:2105.07810 [astro-ph.CO].
  * [40] S. Weinberg, “Quantum contributions to cosmological correlations,” Phys. Rev. D 72 (2005) 043514, arXiv:hep-th/0506236.
  * [41] D. Pirtskhalava, L. Santoni, and E. Trincherini, “Constraints on Single-Field Inflation,” JCAP 1606 (2016) no. 06, 051, arXiv:1511.01817 [hep-th].
  * [42] G. Panagopoulos and E. Silverstein, “Multipoint correlators in multifield cosmology,” arXiv:2003.05883 [hep-th].
  * [43] L. Senatore and M. Zaldarriaga, “On Loops in Inflation,” JHEP 12 (2010) 008, arXiv:0912.2734 [hep-th].
  * [44] S. Melville and E. Pajer, “Cosmological Cutting Rules,” JHEP 05 (2021) 249, arXiv:2103.09832 [hep-th].
  * [45] J. Fumagalli, S. Renaux-Petel, and L. T. Witkowski, “Oscillations in the stochastic gravitational wave background from sharp features and particle production during inflation,” arXiv:2012.02761 [astro-ph.CO].
  * [46] M. Libanov, S. Mironov, and V. Rubakov, “Generalized Galileons: instabilities of bouncing and Genesis cosmologies and modified Genesis,” JCAP 08 (2016) 037, arXiv:1605.05992 [hep-th].
  * [47] T. Kobayashi, “Generic instabilities of nonsingular cosmologies in Horndeski theory: A no-go theorem,” Phys. Rev. D 94 (2016) no. 4, 043511, arXiv:1606.05831 [hep-th].
  * [48] D. Pirtskhalava, L. Santoni, E. Trincherini, and P. Uttayarat, “Inflation from Minkowski Space,” JHEP 12 (2014) 151, arXiv:1410.0882 [hep-th].
  * [49] P. Creminelli, D. Pirtskhalava, L. Santoni, and E. Trincherini, “Stability of Geodesically Complete Cosmologies,” JCAP 11 (2016) 047, arXiv:1610.04207 [hep-th].
  * [50] Y. Cai, Y. Wan, H.-G. Li, T. Qiu, and Y.-S. Piao, “The Effective Field Theory of nonsingular cosmology,” JHEP 01 (2017) 090, arXiv:1610.03400 [gr-qc].
  * [51] J.-W. Chen, J. Liu, H.-L. Xu, and Y.-F. Cai, “Tracing Primordial Black Holes in Nonsingular Bouncing Cosmology,” Phys. Lett. B 769 (2017) 561–568, arXiv:1609.02571 [gr-qc].
  * [52] G. Goon, K. Hinterbichler, A. Joyce, and M. Trodden, “Aspects of Galileon Non-Renormalization,” JHEP 11 (2016) 100, arXiv:1606.02295 [hep-th].
  * [53] A. Nicolis, R. Rattazzi, and E. Trincherini, “The Galileon as a local modification of gravity,” Phys. Rev. D 79 (2009) 064036, arXiv:0811.2197 [hep-th].
  * [54] M. A. Luty, M. Porrati, and R. Rattazzi, “Strong interactions and stability in the DGP model,” JHEP 09 (2003) 029, arXiv:hep-th/0303116.
  * [55] J. Noller, L. Santoni, E. Trincherini, and L. G. Trombetta, “Black Hole Ringdown as a Probe for Dark Energy,” Phys. Rev. D 101 (2020) 084049, arXiv:1911.11671 [gr-qc].
  * [56] J. Gleyzes, D. Langlois, F. Piazza, and F. Vernizzi, “Essential Building Blocks of Dark Energy,” JCAP 08 (2013) 025, arXiv:1304.4840 [hep-th].
  * [57] C. Deffayet, O. Pujolas, I. Sawicki, and A. Vikman, “Imperfect Dark Energy from Kinetic Gravity Braiding,” JCAP 10 (2010) 026, arXiv:1008.0048 [hep-th].
  * [58] N. Arkani-Hamed and J. Maldacena, “Cosmological Collider Physics,” arXiv:1503.08043 [hep-th].
  * [59] E. Pajer, “Building a Boostless Bootstrap for the Bispectrum,” JCAP 01 (2021) 023, arXiv:2010.12818 [hep-th].
  * [60] M. Taoso and A. Urbano, “Non-gaussianities for primordial black hole formation,” arXiv:2102.03610 [astro-ph.CO].
  * [61] S. Garcia-Saenz, S. Renaux-Petel, and J. Ronayne, “Primordial fluctuations and non-Gaussianities in sidetracked inflation,” JCAP 07 (2018) 057, arXiv:1804.11279 [astro-ph.CO].
  * [62] G. A. Palma, S. Sypsas, and C. Zenteno, “Seeding primordial black holes in multifield inflation,” Phys. Rev. Lett. 125 (2020) no. 12, 121301, arXiv:2004.06106 [astro-ph.CO].
  * [63] A. Achucarro, J.-O. Gong, S. Hardeman, G. A. Palma, and S. P. Patil, “Features of heavy physics in the CMB power spectrum,” JCAP 01 (2011) 030, arXiv:1010.3693 [hep-ph].
  * [64] S. Cespedes, V. Atal, and G. A. Palma, “On the importance of heavy fields during inflation,” JCAP 05 (2012) 008, arXiv:1201.4848 [hep-th].
  * [65] A. Achucarro, J.-O. Gong, S. Hardeman, G. A. Palma, and S. P. Patil, “Effective theories of single field inflation when heavy fields matter,” JHEP 05 (2012) 066, arXiv:1201.6342 [hep-th].
  * [66] A. Achucarro, V. Atal, S. Cespedes, J.-O. Gong, G. A. Palma, and S. P. Patil, “Heavy fields, reduced speeds of sound and decoupling during inflation,” Phys. Rev. D 86 (2012) 121301, arXiv:1205.0710 [hep-th].
  * [67] C. Gordon, D. Wands, B. A. Bassett, and R. Maartens, “Adiabatic and entropy perturbations from inflation,” Phys. Rev. D 63 (2000) 023506, arXiv:astro-ph/0009131.
  * [68] S. Groot Nibbelink and B. J. W. van Tent, “Density perturbations arising from multiple field slow roll inflation,” arXiv:hep-ph/0011325.
  * [69] S. Groot Nibbelink and B. J. W. van Tent, “Scalar perturbations during multiple field slow-roll inflation,” Class. Quant. Grav. 19 (2002) 613–640, arXiv:hep-ph/0107272.
  * [70] A. Achucarro, J.-O. Gong, S. Hardeman, G. A. Palma, and S. P. Patil, “Mass hierarchies and non-decoupling in multi-scalar field dynamics,” Phys. Rev. D 84 (2011) 043502, arXiv:1005.3848 [hep-th].
  * [71] R. Gwyn, G. A. Palma, M. Sakellariadou, and S. Sypsas, “Effective field theory of weakly coupled inflationary models,” JCAP 04 (2013) 004, arXiv:1210.3020 [hep-th].


# Computational and Exploratory Landscape Analysis
of the GKLS Generator

Jakub Kudela  0000-0002-4372-2105 Institute of Automation and Computer Science
Brno University of TechnologyBrnoCzech Republic<EMAIL_ADDRESS>and
Martin Juricek  0000-0002-7943-8659 Institute of Automation and Computer
Science
Brno University of TechnologyBrnoCzech Republic<EMAIL_ADDRESS>
(2023)

###### Abstract.

The GKLS generator is one of the most used testbeds for benchmarking global
optimization algorithms. In this paper, we conduct both a computational
analysis and the Exploratory Landscape Analysis (ELA) of the GKLS generator.
We utilize both canonically used and newly generated classes of GKLS-generated
problems and show their use in benchmarking three state-of-the-art methods
(from evolutionary and deterministic communities) in dimensions 5 and 10. We
show that the GKLS generator produces “needle in a haystack” type problems
that become extremely difficult to optimize in higher dimensions. Furthermore,
we conduct the ELA on the GKLS generator and then compare it to the ELA of two
other widely used benchmark sets (BBOB and CEC 2014), and discuss the
meaningfulness of the results.

Benchmarking, Exploratory Landscape Analysis, GKLS, Global optimization,
Black-box optimization

††journalyear: 2023††copyright: rightsretained††ccs: Computing methodologies
Randomized search††ccs: Computing methodologies Heuristic function
construction††ccs: Computing methodologies Continuous space search

##  1\. Introduction

Solving real-world black-box optimization problems is a very challenging task,
even when one has domain knowledge and experience (Long et al., 2022),
especially in problems with expensive function evaluation that require
simulation runs (Kudela and Matousek, 2022). The development, utilization, and
comparison of optimization methods on such real-world problems are usually
prohibitively expensive. For such tasks, benchmarking optimization methods on
artificially constructed testbeds become pivotal, with the expectation that
the behavior of the methods on these benchmark sets translates well into real-
world problems.

Over the years, various benchmark suites have been proposed, in which
different global function properties are represented, such as multi-modality,
separability, ill-conditioning, and various other types of global structures.
In the evolutionary computation community, the two most utilized benchmark
sets are the Black-Box Optimization Benchmarking (BBOB) suite (Hansen et al.,
2012) which is now part of the COCO platform (Hansen et al., 2021), and the
suites that were presented at the Congress on Evolutionary Computation (CEC)
competitions (which started in 2005 and continue to this day) (Kudela, 2022).
As was shown in (Garden and Engelbrecht, 2014), the characteristics of the
functions used in these two benchmarks are quite different. The CEC benchmarks
are constructed by using similar subfunctions, which possibly gives an
advantage to methods that perform well on these fewer subfunctions. It was
also found that the CEC functions share more similarities among themselves
than with those found in the BBOB (Garden and Engelbrecht, 2014).

In the global (deterministic) optimization community, one of the most popular
benchmark sets is the one produced by the GKLS generator (Gaviano et al.,
2003). The GKLS generator constructs classes of test functions (either non-
differentiable, differentiable, or twice-differentiable) for multi-modal,
multi-dimensional box-constrained global optimization. The advantage of the
GKLS generator is that for each generated problem, the location and function
value of its local and global minima are known. Although the GKLS generator
can be used to create various types of problems (based on input parameters),
there are 8 classes of problems (2 for dimensions 2, 3, 4, and 5) each
containing 100 functions that are generally used (Sergeyev and Kvasov, 2006;
Sergeyev et al., 2018; Paulavičius et al., 2020; Stripinis and Paulavičius,
2022). The GKLS generator was also recently used for the construction of
nonlinear model predictive control (Xu and Xu, 2023) and general-constrained
(Sergeyev et al., 2021) test problems.

In order to quantify the low-level properties of optimization problems,
various features of the landscape can be computed (Long et al., 2022). Such
analysis falls under the field of Exploratory Landscape Analysis (ELA)
(Mersmann et al., 2010). The landscape features try to approximate different
aspects of the optimization problem, such as its modality, separability, or
whether or not the problem has plateaus. The most notable uses of ELA are in
the visualization of the problem space of various optimization benchmark
problem sets (Škvorc et al., 2020), and in automated algorithm selection
(Kerschke and Trautmann, 2019a).

In this paper, we conduct computational analysis and ELA of the GKLS-generated
problems. The rest of the paper is structured as follows. Section 2 briefly
describes the GKLS generator. In Section 3, we conduct the computational
analysis of benchmarking three state-of-the-art methods on three classes of
GKLS-generated suits (two canonical and one newly proposed) in dimensions 5
and 10. In Section 4 the ELA of the GKLS-generated suits is presented, along
with the comparison of ELA of the BBOB and CEC 2014 suits. Finally,
conclusions are drawn in Section 5.

|  |
---|---|---
$h=2$, type D, $d=0.9$, $r=0.1$ |  $h=5$, type ND, $d=0.66$, $r=0.5$ |  $h=10$, type D, $d=0.66$, $r=0.3$
|  |
$h=20$, type D, $d=0.9$, $r=0.1$ |  $h=50$, type D, $d=0.66$, $r=0.2$ |  $h=100$, type ND, $d=0.66$, $r=0.5$
Figure 1. Functions generated by the GKLS in dimension $D=2$ ($f^{*}=-1$ for
all functions).

##  2\. GKLS generator

In the GKLS generator, a prespecified number of test problems (a class of
problems) is constructed by defining a convex quadratic function (a
paraboloid) which is systematically distorted by polynomials in order to
produce local (and one global) minima. The input parameters for this
construction are the following: type of the problem (ND: non-differentiable,
D: differentiable, D2: twice-differentiable) problem dimension ($D$), number
of local minima ($h$), the value of the global minimum ($f^{*}$), radius ($r$)
of the attraction region of the global minimizer, and the distance ($d$) from
the global minimizer to the vertex of the quadratic function. All problems are
constructed on $[-1,1]^{D}$.

A visualization of different functions that can be generated by the GKLS and
the effect of different parameter choices is shown in Figure 1. Several
interesting observations can be made regarding the generated functions.
Firstly, they are all relatively well-conditioned. The local minimum of the
“big” paraboloid is always in the domain and has a function value of 0. The
“attraction regions” of the different local minima do not overlap (this is by
design) - this also means that they become more shallow when their number
increases.

The eight most used classes (each with 100 functions) of GKLS-generated
problems are shown in Table 1. An interesting thing to note is that these
classes are quite similar, with type, $f^{*}$, and $h$ parameters being the
same. Another feature of the GKLS generator is that it has a built-in
pseudorandom number generator. If the same parameters are given as an input to
GKLS an identical class of functions will be produced. However, the
pseudorandom number generator (or, more precisely, its seed) does not depend
on all input parameters but is mainly dependent on $h$. This has some
interesting effects. For instance, two problems in different classes in the
same dimension and with the same number (i.e., problem number 1 in class 1 and
problem number 1 in class 2), and with the same value of $h$ will have the
same location of the local minimum of the “big” paraboloid. If the classes
also share the same $d$, the problems with the same number will also have the
same location of the global minimum. This issue can be seen in problem classes
1 and 2, and 7 and 8.

Table 1. Most used GKLS test classes. Class | “Difficulty” | Type | $D$ | $f^{*}$ | $d$ | $r$ | $h$
---|---|---|---|---|---|---|---
1 | simple | D | 2 | -1 | 0.90 | 0.20 | 10
2 | hard | D | 2 | -1 | 0.90 | 0.10 | 10
3 | simple | D | 3 | -1 | 0.66 | 0.20 | 10
4 | hard | D | 3 | -1 | 0.90 | 0.20 | 10
5 | simple | D | 4 | -1 | 0.66 | 0.20 | 10
6 | hard | D | 4 | -1 | 0.90 | 0.20 | 10
7 | simple | D | 5 | -1 | 0.66 | 0.30 | 10
8 | hard | D | 5 | -1 | 0.66 | 0.20 | 10
|  |
---|---|---
AGSK | LSHADE | BIRMIN
Figure 2. Convergence plots of the three methods in dimension $D=5$. |  |
---|---|---
GKLS simple | GKLS hard | GKLS mod
Figure 3. Empirical cumulative distribution of simulated runtimes, measured in
number of $f$-evaluations for the $51$ targets $10^{[-8..2]}$ in dimension
$D=5$.

##  3\. Computational Analysis

Although the classes in the same dimension might share some characteristics,
what is arguably more important is whether it means that different algorithms
will “perform” in the same way on these classes of problems. To study this
effect, we set up the following computational analysis. We run three state-of-
the-art methods both from the evolutionary computation (EC) and deterministic
optimization communities on the “canonical” GKLS-generated problems in
dimensions 5 and 10 (which uses the same “simple” and “hard” parameters as
dimension 5). We also construct a new class (“mod”) 50 of GKLS-generated
problems (again, in $D=5$, and $10$) by the following procedure:

  * •

Each problem $i=1,\dots,50$ in this class will have different values of the
parameters $\text{type}_{i}$, $d_{i}$, $r_{i}$, and $h_{i}$, but the same
$f_{i}^{*}=-1$.

  * •

$\text{type}_{i}$ is decided by a coin flip between types D and ND (with the
same probability).

  * •

$d_{i}$ is a uniformly distributed random number on [0,1].

  * •

$r_{i}=d_{i}/u_{i}$, where $u_{i}$ is a uniformly distributed random integer
on [2,10], i.e. $r_{i}\in[d_{i}/2,d_{i}/10]$.

  * •

$h_{i}=\text{round}(10^{c_{i}})$, where $c_{i}$ is a uniformly distributed
random number on [1,3], i.e. $h\in[10,10^{3}]$.

From the EC side, we chose two methods to run on the GKLS-generated problems.
The first selected method was Adaptive Gaining-Sharing Knowledge (AGSK), which
was the runner-up of the CEC’20 competition. The algorithm improved the
original GSK algorithm by adding adaptive settings to its two control
parameters: the knowledge factor and ratio, which control junior and senior
gaining and sharing phases during the optimization process (Mohamed et al.,
2020). The second method is L-SHADE or Success-history based adaptive
differential evolution with linear population size reduction (Tanabe and
Fukunaga, 2014). This metaheuristic method has its basis in adaptive
differential evolution, which involves success-history-based parameter
adaptation. The proposed method then provides an extension in the form of
using linear population size reduction, which results in population size
reduction according to a linear function. From the deterministic methods, we
selected the best-performing method from a recent extensive numerical study
(Stripinis and Paulavičius, 2022), which evaluated 64 derivative-free
algorithms on the test problems from the DIRECTGOLib (Stripinis and
Paulavičius, 2022) and on the GKLS generator. The chosen method was BIRMIN
(Paulavičius et al., 2020), which is a globally-biased hybridized version of
the BIRECT method (Paulavičius et al., 2018), which was also “trained” on the
GKLS generator.

|  |
---|---|---
AGSK | LSHADE | BIRMIN
Figure 4. Convergence plots of the three methods in dimension $D=10$. |  |
---|---|---
GKLS simple | GKLS hard | GKLS mod
Figure 5. Empirical cumulative distribution of simulated runtimes, measured in
number of $f$-evaluations for the $51$ targets $10^{[-8..2]}$ in dimension
$D=10$.

The implementation and parameter choices for LSHADE and AGSK were taken from
(Biswas et al., 2021) (the implementations are available at the
GitHub111https://github.com/subhodipbiswas/MadDE of one of the authors). The
implementation and parameter choices for BIRMIN were taken from the
DIRECTGOLib (Stripinis and Paulavičius, 2022).

For the numerical comparison, we chose to run each of the three methods once
on every problem from each of the three classes (100 problems in both “simple”
and “hard” classes, and 50 problems in the “mod” class) in dimensions
$D=[5,10]$, with a budget of $5\cdot 10^{4}\cdot D$ available function
evaluations. For every run, if the objective function value of the resulting
solution was less than or equal to 1E–8, it was considered as zero. All
algorithms were run in a MATLAB R2022a, on a PC with a 3.2 GHz Core I5 CPU, 16
GB RAM, and Windows 10. The code for the experiments (and the generator for
the test problems) can be found at the author’s
Github222https://github.com/JakubKudela89/GKLS-GECCO.

We begin the discussion on the results of the computations on classes in
dimension $D=5$. The convergence plots of the three methods on the three
classes are shown in Figure 2, while the Empirical cumulative distributions
(ECDs) of simulated runtimes, measured in the number of function evaluations
for 51 targets $10^{[-8..2]}$ (similar analysis which is done in the COCO
platform) are shown in Figure 3. The first thing one can find in these results
is that although the “simple” and “hard” classes are composed of similar
problems (as discussed above), there is a noticeable difference in the
behavior of the three algorithms. On the “simple” class all three algorithms
were able to find either a good or the optimal solutions faster than on the
“hard” class. There is also a quite large difference between the performance
of the three different methods - BIRMIN clearly dominated the two EC methods
(and was able to find the optimal solution for a large portion of the
problems), and AGSK turned out to be better at finding good solutions at the
later stages of the search than LSHADE. However, both AGSK and LSHADE
experienced a plateau in their convergence between $10^{4}$ and $4\cdot
10^{4}$ function evaluations and got stuck in local optima. These methods
would undoubtedly benefit from using restart strategies (Tanabe and Fukunaga,
2015) or from hybridization with model-based optimization (Okulewicz and
Zaborski, 2021), but we did not pursue these possibilities.

The results change quite dramatically when looking at the “mod” class.
Although the performance of BIRMIN is still superior to that of the two EC
methods, the margin narrowed substantially. What is more, the relative
performance (against the “hard” class) of AGSK decreased, while for LSHADE it
increased. For all three classes, the local minimum of the “big” paraboloid
(error value 1) was found by every method within a few hundred function
evaluations for BIRMIN and a few thousand function evaluations for the EC
methods.

|  |
---|---|---
| dimension $D=5$ |
|  |
| dimension $D=10$ |
Figure 6. Relashionship between $r,d,h$, type, and the best function value
found on the “mod” class problems for the three methods for dimensions $D=5$
and $D=10$.

The results for dimension $D=10$, shown in Figures 4 and 5 show another
substantial change. Looking at the ECD plot (notice the change in the range of
the y-axis compared to a similar figure in dimension $D=5$) in Figure 5, we
can find that all three classes are basically equivalent. The plateaus in the
ECD plots between the 0.19 and 0.20 values on the y-axis indicate the points
where the methods found the local minimum of the “big” paraboloid. Although
this happened relatively early for all methods, finding better local optima
proved to be extremely challenging. And the global optimum was found only
twice in the “simple” class (once by BIRMIN and once by AGSK on different
problems) and once in the “hard” class (by AGSK). In the “mod” class, no
method was able to find the global optimum for any of the problems, but there
were more “good” local optima being found, especially by DIRMIN (hence the
slightly better ECD plot).

As the problems “mod” class were generated using different values of the type,
$r,d$, and $h$ parameters, we can study their effect on the solutions obtained
by the three methods. This relationship is depicted in Figure 6. We focus
first on dimension $D=5$. First, we can see that the type of the problem
(differentiable or non-differentiable) had practically no effect on the
quality of the obtained solution. Rather intuitively, larger values of $r$
(region of attraction of the global minimizer) resulted in the methods finding
the global optimum more easily. What is not so intuitive is the dependence on
$d$ (distance from the global minimizer to the local minimum of the “big”
paraboloid). Here, we can see that larger values of $d$ meant that the methods
had a better chance of finding the global optimum. This was probably caused by
the way $r$ was generated (as a fraction of $d$), which meant that problems
with larger values of $d$ had also larger values of $r$. The last parameter we
look at is $h$ (the number of local minima). While the difficulty of finding
the global optimum seems to be unaffected by $h$, the parameter did have an
effect if the best-found solution was not the global optimum. Increasing $h$
had the effect of worsening the best-found solution.

|
---|---
dimension $D=5$ | dimension $D=10$
Figure 7. Relative frequencies of function values of local minima of the different classes. |
---|---
a) | b)
Figure 8. Scatter plots of the dependence of a) function values of local
minima with function value ¡ 0, b) number of local minima with function value
¡ 0 on $h$ and dimension $D$.

We can trace this effect back to the function plots shown in Figure 1, where
having more local minima meant that they became shallower. We can see this
effect in Figure 7, where the relative frequencies of function values of the
local minima for the three classes are plotted (the “simple” and “hard”
classes have almost the same values of all the local minima - this is again, a
consequence of the implemented pseudorandom number generator).

In dimension $D=10$, the dependence of the best-found function value on the
different values of the parameters vanishes. This is in line with the
observation that the behavior of the three methods on three different classes
in dimension $D=10$ almost did not change (Figure 5). However, one would still
expect a similar dependence on $h$ as was seen in dimension $D=5$. Especially
since Figure 7 shows a noticeable shift of the function values of the local
minima to higher values for dimension $D=10$ (when compared to dimension
$D=5$). This is explained in Figure 8, where we can see that although there
are really fewer local minima in dimension $D=10$ with negative function
values, they are spread out more evenly than in the case of dimension $D=5$,
even when the number of local minima is high.

The functions that the GKLS generator produces are of the “needle in a
haystack” kind (Droste et al., 2002). If one does not stumble upon the region
of the space where the global minimizer (or at least a local minimizer with
“good” function value) has its region of attraction, either by chance (in the
case of the EC methods) or by space partition (whose cost is bound to be
exponential in the dimension $D$), the solution one gets is the local minimum
of the “big” paraboloid. The functions give no hints on where such good points
might be.

##  4\. Exploratory Landscape Analysis

We use ELA features to show how the GKLS-generated problems compare to the
BBOB and CEC 2014 benchmark suits. In order to calculate the ELA features, we
used the flacco library (Kerschke and Trautmann, 2019b). We chose ELA feature
sets which only require samples of input and function value pairs: ela_distr,
ela_meta, disp, nbc, pca, and ic, and dimension $D=10$ for all considered
suits. It was recently shown that the ELA features are sensitive to sampling
strategy (Renau et al., 2020) and function transformations (Škvorc et al.,
2022, 2021). We chose to ignore the features that were sensitive to function
transformations and used uniform sampling with $250\cdot D$ samples for the
computation of the ELA features in all suits. There were 24 problems in the
BBOB suit, 30 problems in the CEC 2014 suit, 100 problems in the GKLS “simple”
and “hard” classes, and 50 problems in the GKLS “mod” class.

We then followed the methodology described in (Škvorc et al., 2020) for the
selection and visualization of the relevant ELA features. The features that
produced constant results on every problem and those that produced invalid
values were removed. Another batch of features that got removed were the
highly correlated ones. The 14 features that remained, along with their
maximum and minimum values on the different suits are shown in Table 2. We can
see that the ranges of values of the chosen features for the three GKLS-
generated suits are very similar. The ranges of the feature values of the
GKLS-generated suits are also generally narrower than that of both BBOB and
CEC 2014 suits. Both of these observations should be expected, as (as we have
seen in the previous section) the GKLS generator produces somehow limited
types of functions. One interesting observation can be made regarding the
value of the ela_distr.number_of_peaks on the GKLS “mod” suit. Although in
this suit there are problems with hundreds of local minima, they are
practically all too shallow to be noticeable in the sampled points.

Table 2. Minimum and Maximum values of the relevant ELA features on the
different benchmark sets. Extremal values are highlighted in bold.

| BBOB | CEC 2014 | GKLS simple | GKLS hard | GKLS mod
---|---|---|---|---|---
| min | max | min | max | min | max | min | max | min | max
ela_distr.skewness | -2.97E+00 | 8.28E+00 | -6.63E-01 | 6.47E+00 | 1.02E-01 | 4.18E-01 | 1.89E-01 | 4.22E-01 | 2.01E-01 | 4.46E-01
ela_distr.kurtosis | -4.94E-01 | 9.67E+01 | -3.38E-01 | 6.50E+01 | -3.53E-01 | 2.15E-01 | -3.40E-01 | 1.66E-01 | -3.56E-01 | 2.47E-01
ela_distr.number_of_peaks | 1.00E+00 | 1.80E+01 | 1.00E+00 | 2.60E+01 | 1.00E+00 | 2.00E+00 | 1.00E+00 | 2.00E+00 | 1.00E+00 | 2.00E+00
ela_meta.lin_simple.intercept | -9.17E+02 | 9.62E+08 | 5.22E+02 | 5.63E+10 | 4.88E+00 | 9.00E+00 | 4.80E+00 | 8.92E+00 | 5.06E+00 | 9.00E+00
ela_meta.lin_w_interact.adj_r2 | 2.14E-04 | 1.00E+00 | -9.41E-04 | 9.04E-01 | 6.91E-01 | 8.98E-01 | 6.82E-01 | 8.93E-01 | 7.04E-01 | 8.95E-01
ela_meta.quad_simple.adj_r2 | 3.98E-03 | 1.00E+00 | -3.61E-03 | 9.88E-01 | 9.93E-01 | 1.00E+00 | 9.94E-01 | 1.00E+00 | 9.90E-01 | 1.00E+00
ela_meta.quad_w_interact.adj_r2 | 3.67E-05 | 1.00E+00 | -1.26E-02 | 1.00E+00 | 9.94E-01 | 1.00E+00 | 9.95E-01 | 1.00E+00 | 9.91E-01 | 1.00E+00
disp.ratio_mean_25 | 8.43E-01 | 1.00E+00 | 8.60E-01 | 1.01E+00 | 8.44E-01 | 8.84E-01 | 8.37E-01 | 8.84E-01 | 8.43E-01 | 8.74E-01
disp.ratio_median_02 | 6.48E-01 | 1.06E+00 | 6.57E-01 | 1.02E+00 | 6.21E-01 | 7.22E-01 | 6.25E-01 | 7.19E-01 | 6.32E-01 | 7.20E-01
nbc.nb_fitness.cor | -6.41E-01 | -1.78E-01 | -6.30E-01 | -1.90E-01 | -4.08E-01 | -3.52E-01 | -4.20E-01 | -3.58E-01 | -4.12E-01 | -3.65E-01
pca.expl_var.cov_init | 9.09E-02 | 9.09E-01 | 9.09E-02 | 9.09E-01 | 6.36E-01 | 8.18E-01 | 6.36E-01 | 8.18E-01 | 6.36E-01 | 8.18E-01
pca.expl_var.cor_init | 8.18E-01 | 9.09E-01 | 8.18E-01 | 9.09E-01 | 8.18E-01 | 8.18E-01 | 8.18E-01 | 8.18E-01 | 8.18E-01 | 8.18E-01
pca.expl_var_PC1.cov_init | 1.07E-01 | 1.00E+00 | 1.10E-01 | 1.00E+00 | 4.97E-01 | 7.48E-01 | 4.98E-01 | 7.41E-01 | 5.23E-01 | 7.52E-01
pca.expl_var_PC1.cor_init | 1.02E-01 | 1.83E-01 | 9.94E-02 | 1.76E-01 | 1.66E-01 | 1.85E-01 | 1.66E-01 | 1.82E-01 | 1.64E-01 | 1.80E-01

For further analysis, the values of the ELA features on the different
benchmark sets were normalized, and we used Principal Component Analysis (PCA)
to reduce the number of features even further. Figure 9 shows a representation
of the 14 PCA components obtained when comparing the ELA features (normalized)
calculated on the combined set of CEC 2014, BBOB, and all GKLS-generated
problems. Using the first 7 components explained 99.68% of the variance. For
visualizing the results, we used the t-Distributed Stochastic Neighbor
Embedding (t-sne). In this visualization, which is shown in Figure 10,
benchmark problems that have similar ELA features should be shown close to
each other.

We can see that the t-sne visualization grouped most of the functions from the
BBOB and CEC 2014 suits together (in a few groups), while the “similar”
problems generated in the three GKLS suites take up most of the space. It
should be pointed out that the t-sne procedure is random, producing different
plots every time. However, in our experiments, all plots looked qualitatively
the same as Figure 10. Interestingly, some of the GKLS-generated problems show
common traits to a few problems in the other two benchmark suits. The three
BBOB functions that are close to the GKLS-generated ones ($f_{3}$, $f_{5}$,
and $f_{13}$) are shown (in dimension $D=2$) in Figure 11. All three are
probably close to some instance of the GKLS-generated “big” paraboloid, which
would explain the closeness of the problems in the t-sne plot.

Figure 9. The amount of explained variance per component when performing PCA on the ELA features calculated on the combined set of all benchmark problems. Figure 10. The t-sne visualization of the ELA features (after normalization and using the first seven components from the PCA) of the benchmark sets. |  |
---|---|---
$f_{3}$ (Rastrigin separable) |  $f_{5}$ (Linear slope) |  $f_{13}$ (Sharp ridge)
Figure 11. The three functions from the BBOB suit that are the “closest” to
some GKLS-generated instances.

##  5\. Conclusion

In this paper, we analyzed the problems constructed by the GKLS generator. In
the computational analysis, we have shown that it produces what are basically
“needle in a haystack ” problems which get extremely difficult to optimize as
the problem dimension grows. The GKLS generator could be successfully used for
benchmarking state-of-the-art methods in lower dimensions ($D=5$) on some of
the simpler instances. However, in the higher dimension ($D=10$), the
performance of the three considered methods was hard to differentiate as the
problems became extremely difficult for the given computational budget. This
difficulty of the generated instances was also largely unaffected by the
choice of parameters that the generator has. Although increasing the
computational budget might bring additional insight, different restart
strategies for the EC methods would have to be used, as they both plateaued
long before reaching the maximum available function evaluations.

The type of function generated by GKLS (differentiable or non-differentiable)
had practically no effect on its “usefulness” in benchmarking. It is possible
that the GKLS generator could be modified to have a much “deeper” local
minima. As the task of finding the global minimum is practically impossible in
higher dimensions, having problems with lots of “good” local minima (i.e.,
better ones than the local minimum of the “big” paraboloid) could be useful
for analyzing the exploration capabilities of optimization methods.

It is not very clear how one could meaningfully use the results of the
computations of the ELA features or the t-sne plots on such “needle in a
haystack” problems. The superficial closeness the some of the GKLS-generated
problems to the three BBOB problems, and the large variation in the t-sne
dimensions of the GKLS-generated problems hint at the limits where ELA can be
meaningfully used. It is probably impossible to have any sample-based features
that would both uncover that the problem is a “needle in a haystack” and be
computationally tractable (as it would amount to finding the “needle” in a
reasonable amount of function evaluations).

Lastly, it would be interesting to find real-world black-box continuous
problems which have the “needle in a haystack” character and show if the
methods that were “trained” on such problems (such as BIRMIN) really offer an
advantage (at least in lower dimensions).

###### Acknowledgements.

This work was supported by the IGA BUT No. FSI-S-23-8394 “Artificial
intelligence methods in engineering tasks”.

## References

  * (1)
  * Biswas et al. (2021) Subhodip Biswas, Debanjan Saha, Shuvodeep De, Adam D Cobb, Swagatam Das, and Brian A Jalaian. 2021\.  Improving differential evolution through Bayesian hyperparameter optimization. In _2021 IEEE Congress on Evolutionary Computation (CEC)_. IEEE, 832–840.
  * Droste et al. (2002) Stefan Droste, Thomas Jansen, and Ingo Wegener. 2002\.  Optimization with randomized search heuristics—the (A) NFL theorem, realistic scenarios, and difficult functions.  _Theoretical Computer Science_ 287, 1 (2002), 131–144.
  * Garden and Engelbrecht (2014) Robert W Garden and Andries P Engelbrecht. 2014.  Analysis and classification of optimisation benchmark functions and benchmark suites. In _2014 IEEE Congress on Evolutionary Computation (CEC)_. IEEE, 1641–1649.
  * Gaviano et al. (2003) Marco Gaviano, Dmitri E Kvasov, Daniela Lera, and Yaroslav D Sergeyev. 2003.  Algorithm 829: Software for generation of classes of test functions with known local and global minima for global optimization.  _ACM Transactions on Mathematical Software (TOMS)_ 29, 4 (2003), 469–480.
  * Hansen et al. (2012) Nikolaus Hansen, Anne Auger, Steffen Finck, and Raymond Ros. 2012.  _Real-Parameter Black-Box Optimization Benchmarking 2012: Experimental Setup_.  Technical Report. INRIA.
  * Hansen et al. (2021) N. Hansen, A. Auger, R. Ros, O. Mersmann, T. Tušar, and D. Brockhoff. 2021\.  COCO: A Platform for Comparing Continuous Optimizers in a Black-Box Setting.  _Optimization Methods and Software_ 36 (2021), 114–144.  Issue 1. https://doi.org/10.1080/10556788.2020.1808977
  * Kerschke and Trautmann (2019a) Pascal Kerschke and Heike Trautmann. 2019a.  Automated algorithm selection on continuous black-box problems by combining exploratory landscape analysis and machine learning.  _Evolutionary computation_ 27, 1 (2019), 99–127.
  * Kerschke and Trautmann (2019b) Pascal Kerschke and Heike Trautmann. 2019b.  Comprehensive feature-based landscape analysis of continuous and constrained optimization problems using the R-package flacco.  In _Applications in Statistical Computing_. Springer, 93–123.
  * Kudela (2022) Jakub Kudela. 2022\.  A critical problem in benchmarking and analysis of evolutionary computation methods.  _Nature Machine Intelligence_ 4 (2022), 1238–1245.
  * Kudela and Matousek (2022) Jakub Kudela and Radomil Matousek. 2022.  Recent advances and applications of surrogate models for finite element method computations: a review.  _Soft Computing_ 26, 24 (2022), 13709–13733.
  * Long et al. (2022) Fu Xing Long, Diederick Vermetten, Bas van Stein, and Anna V Kononova. 2022.  BBOB Instance Analysis: Landscape Properties and Algorithm Performance across Problem Instances.  _arXiv preprint arXiv:2211.16318_ (2022).
  * Mersmann et al. (2010) Olaf Mersmann, Mike Preuss, and Heike Trautmann. 2010\.  Benchmarking Evolutionary Algorithms: Towards Exploratory Landscape Analysis. In _Parallel Problem Solving from Nature, PPSN XI: 11th International Conference, Kraków, Poland, September 11-15, 2010, Proceedings, Part I 11_. Springer, 73–82.
  * Mohamed et al. (2020) Ali Wagdy Mohamed, Anas A Hadi, Ali Khater Mohamed, and Noor H Awad. 2020.  Evaluating the performance of adaptive gainingsharing knowledge based algorithm on CEC 2020 benchmark problems. In _2020 IEEE Congress on Evolutionary Computation (CEC)_. IEEE, 1–8.
  * Okulewicz and Zaborski (2021) Michał Okulewicz and Mateusz Zaborski. 2021.  Benchmarking SHADE algorithm enhanced with model based optimization on the BBOB noiseless testbed. In _Proceedings of the Genetic and Evolutionary Computation Conference Companion_. 1259–1266.
  * Paulavičius et al. (2018) Remigijus Paulavičius, Lakhdar Chiter, and Julius Žilinskas. 2018.  Global optimization based on bisection of rectangles, function values at diagonals, and a set of Lipschitz constants.  _Journal of Global Optimization_ 71 (2018), 5–20.
  * Paulavičius et al. (2020) Remigijus Paulavičius, Yaroslav D Sergeyev, Dmitri E Kvasov, and Julius Žilinskas. 2020\.  Globally-biased BIRECT algorithm with local accelerators for expensive global optimization.  _Expert Systems with Applications_ 144 (2020), 113052.
  * Renau et al. (2020) Quentin Renau, Carola Doerr, Johann Dreo, and Benjamin Doerr. 2020.  Exploratory landscape analysis is strongly sensitive to the sampling strategy. In _Parallel Problem Solving from Nature–PPSN XVI: 16th International Conference, PPSN 2020, Leiden, The Netherlands, September 5-9, 2020, Proceedings, Part II 16_. Springer, 139–153.
  * Sergeyev et al. (2018) Ya D Sergeyev, DE Kvasov, and MS Mukhametzhanov. 2018.  On the efficiency of nature-inspired metaheuristics in expensive global optimization with limited budget.  _Scientific reports_ 8, 1 (2018), 1–9.
  * Sergeyev and Kvasov (2006) Yaroslav D Sergeyev and Dmitri E Kvasov. 2006.  Global search based on efficient diagonal partitions and a set of Lipschitz constants.  _SIAM Journal on Optimization_ 16, 3 (2006), 910–937.
  * Sergeyev et al. (2021) Yaroslav D Sergeyev, Dmitri E Kvasov, and Marat S Mukhametzhanov. 2021.  A Generator of Multiextremal Test Classes With Known Solutions for Black-Box-Constrained Global Optimization.  _IEEE Transactions on Evolutionary Computation_ 26, 6 (2021), 1261–1270.
  * Škvorc et al. (2020) Urban Škvorc, Tome Eftimov, and Peter Korošec. 2020\.  Understanding the problem space in single-objective numerical optimization using exploratory landscape analysis.  _Applied Soft Computing_ 90 (2020), 106138.
  * Škvorc et al. (2021) Urban Škvorc, Tome Eftimov, and Peter Korošec. 2021\.  The effect of sampling methods on the invariance to function transformations when using exploratory landscape analysis. In _2021 IEEE Congress on Evolutionary Computation (CEC)_. IEEE, 1139–1146.
  * Škvorc et al. (2022) Urban Škvorc, Tome Eftimov, and Peter Korošec. 2022\.  A Comprehensive Analysis of the Invariance of Exploratory Landscape Analysis Features to Function Transformations. In _2022 IEEE Congress on Evolutionary Computation (CEC)_. IEEE, 1–8.
  * Stripinis and Paulavičius (2022) Linas Stripinis and Remigijus Paulavičius. 2022\.  An extensive numerical benchmark study of deterministic vs. stochastic derivative-free global optimization algorithms.  _arXiv preprint arXiv:2209.05759_ (2022).
  * Stripinis and Paulavičius (2022) Linas Stripinis and Remigijus Paulavičius. 2022.  _DIRECTGOLib - Global Optimization test problems Library_.  https://doi.org/10.5281/zenodo.6617799
  * Tanabe and Fukunaga (2015) Ryoji Tanabe and Alex Fukunaga. 2015.  Tuning differential evolution for cheap, medium, and expensive computational budgets. In _2015 IEEE Congress on Evolutionary Computation (CEC)_. IEEE, 2018–2025.
  * Tanabe and Fukunaga (2014) Ryoji Tanabe and Alex S Fukunaga. 2014.  Improving the search performance of SHADE using linear population size reduction. In _2014 IEEE congress on evolutionary computation (CEC)_. IEEE, 1658–1665.
  * Xu and Xu (2023) JiaHong Xu and LiHong Xu. 2023.  N-GKLS: An NMPC problem generator and a test platform for NMPC solvers.  _Expert Systems with Applications_ 214 (2023), 119029.


# Toward Super-Resolution for Appearance-Based Gaze Estimation

Galen O’Shea
Carleton Unviersity
Ottawa, ON, Canada
<EMAIL_ADDRESS>Majid Komeili
Carleton Unviersity
Ottawa, ON, Canada
<EMAIL_ADDRESS>
###### Abstract

Gaze tracking is a valuable tool with a broad range of applications in various
fields, including medicine, psychology, virtual reality, marketing, and
safety. Therefore, it is essential to have gaze tracking software that is
cost-efficient and high-performing. Accurately predicting gaze remains a
difficult task, particularly in real-world situations where images are
affected by motion blur, video compression, and noise. Super-resolution has
been shown to improve image quality from a visual perspective. This work
examines the usefulness of super-resolution for improving appearance-based
gaze tracking. We show that not all SR models preserve the gaze direction. We
propose a two-step framework based on SwinIR super-resolution model. The
proposed method consistently outperforms the state-of-the-art, particularly in
scenarios involving low-resolution or degraded images. Furthermore, we examine
the use of super-resolution through the lens of self-supervised learning for
gaze prediction. Self-supervised learning aims to learn from unlabelled data
to reduce the amount of required labeled data for downstream tasks. We propose
a novel architecture called “SuperVision” by fusing an SR backbone network to
a ResNet18 (with some skip connections). The proposed SuperVision method uses
5x less labeled data and yet outperforms, by 15%, the state-of-the-art method
of GazeTR which uses 100% of training data. We will make our code publicly
available upon publication.

##  1 Introduction

Gaze estimation has been regarded as an valuable research area with a broad
range of applications. Research into gaze tracking has often been divided into
model-based and appearance-based approaches. While model-based methods rely on
specialized hardware, appearance-based methods only require a camera [17].
Although there have been improvements in appearance-based gaze tracking,
accurately predicting gaze remains a difficult task, particularly outside
controlled environments and curated datasets [1, 48, 39, 14]. Furthermore,
research in appearance-based gaze tracking compensate for the lack of
specialized hardware by using large and complex machine learning algorithms to
extract information [51, 2, 3], leaving the data overlooked.

Figure 1: Top: an example of super-resolution improving gaze prediction.
Bottom: an example of super-resolution degrading the gaze prediction by
hallucinating a frontal gaze. The blue line is the gaze prediction and and the
red line is the ground truth. Super-resolution is done using SwinIR [24] and
GFP-GAN [43] for the top and the bottom rows respectively. Gaze Estimation is
done using Full-Face [52].

Super-resolution is a technique used to increase the resolution of an image
beyond the original image resolution. Super-resolution has been shown to
improve image quality from a visual perspective [23, 45, 44, 43, 24, 35]. For
instance, SR can increase detail and clarity while removing degradations such
as motion blur, file compression, and pixelation. Super-resolution has been
shown to be useful in many applications, including medical imaging [36],
remote sensing [25], surveillance and video processing [26]. While SR has been
shown to create visually appealing images, its usefulness for gaze prediction
has not been investigated. Figure 1 shows SR as a preprocessing step for gaze
prediction. We rigorously examine the usefulness of SR for gaze prediction.

Furthermore, we examine the use of SR through the lens of self-supervised
learning for gaze prediction. Self-supervised learning allows a machine
learning model to learn from unlabelled data which can be easier and less
expensive to acquire than labeled data. For example, a self-supervised model
might be trained to predict the rotation or color of an image which can be
done without the need for manual annotation or labelling. Once the self-
supervised model is trained, it can be fine-tuned for a specific task using a
smaller amount of labeled data. Gaze prediction, like many traditional
supervised learning tasks, requires a large amount of labeled data. A sample
efficient model opens up possibility of building gaze models for scenarios
where obtaining a large labelled gaze dataset is challenging, for example gaze
tracking in infants [13], older adults [4] or various types of animals. Gaze
tracking in animals can be useful for studying their behaviour, cognition, and
perception [47, 33, 16].

Contributions of this paper are as follows.

  1. 1.

We rigorously examine the usefulness of SR when used as a preprocessing step
for appearance-based gaze prediction. We show that not all SR models preserve
the gaze direction.

  2. 2.

We propose a two-step framework based on SR and achieve state-of-the-art
results on the MPIIFaceGaze [51] dataset. Additionally, the performance of the
proposed method was evaluated on low-resolution images and various
degradations, and the proposed method consistently outperformed other methods.

  3. 3.

We propose a novel architecture called “SuperVision” by fusing an SR backbone
network to a ResNet (with some skip connections). The proposed SuperVision
method uses only 20% of the labeled data and yet outperforms, by 15%, the
state-of-the-art method of GazeTR [10] which uses 100% of data.

##  2 Related Work

###  2.1 Traditional Approach

In the field of gaze tracking research, there are two primary approaches:
model-based and appearance-based methods [17]. Model-based methods aim to
predict gaze by using either a geometric 3D eye model that analyzes the
corneal reflection [32], or shape-based methods that utilize the pupil centre
[42] or contours of the iris [14]. Advancements in corneal-reflection
techniques have enabled the transition from fixed head positions [31] to
multiple head poses and lighting conditions [54]. However, due to the
complexity of these methods, researchers often employ costly and specialized
hardware (e.g., depth sensors) that may not be practical for general-purpose
gaze tracking [1, 48, 39, 14]. Additionally, although these methods exhibit
outstanding performance in a controlled laboratory setting, they are less
reliable in low light conditions and unconstrained environments.

Appearance-based methods frame gaze estimation as a regression problem by
mapping gaze images to a corresponding gaze vector, while avoiding the need
for specialized hardware and using only a camera. However, early appearance-
based methods required time-consuming head-pose calibration for each
participant, which led to research on reducing the number of training examples
using semi-supervised Gaussian regression methods [46] and finding an optimal
set of training samples using adaptive linear regression [28]. Despite these
efforts, calibration remained insufficient as models failed to generalize to
new subjects and head poses, leading to research on solving head pose and
subject-related issues using a pose-based clustering method [38] and
compensating for bias via regression [27]. To address free head motion, eye
image synthesis was later employed [29]. Generalization problems were further
handled using cross-subject training methods [30] and learning-by-synthesis
methods [37].

###  2.2 Deep Learning Approach

Early approaches to gaze tracking faced significant challenges in adapting to
new subjects and positions despite efforts to enhance their performance. As a
result, research shifted to deep learning approaches such as the Multimodal
Convolutional Neural Network (CNN), which concatenate eye images with a head
pose estimate [51]. In subsequent work, the authors designed a VGG-inspired
architecture that extended their previous work [53]. Other methods included
dilated-convolutions to extract high-level features without reducing spatial
resolution [6]. In related studies on multi-stream CNNs, researchers adopted
data fusion, which involved merging datasets while maintaining separate
validation sets to evaluate the performance of each dataset independently [2].
Although state-of-the-art performance was achieved, the variation in accuracy
across datasets was attributed to differences in the resolution of the eye
patch. The researchers suggested that deep learning models for gaze estimation
could benefit from higher accuracy if the datasets had more pixels within each
eye patch [2].

####  2.2.1 Attention Approach

Studies have indicated that combining eye feature maps can enhance the
accuracy of gaze tracking. For instance, a study utilized a VGG-16 network per
eye, with concatenated downstream features [11]. Other works, such as [6],
employed quad-stream architecture to extract singular and joint features from
both eyes. Subsequently, researchers explored the use of the attention
mechanism to extract joint eye features. In one novel approach, [3] proposed
an Adaptive Feature Fusion technique that stacked eye feature maps based on
their similarity, using a self-attention mechanism. More recently, transformer
architectures have emerged as dominant in the field, including for gaze
tracking. The current state-of-the-art model, GazeTR, uses a hybrid CNN and
transformer architecture [9].

####  2.2.2 Generative Adversarial Networks Approach

In the past, the predominant approach to solving the gaze estimation problem
has been to develop complex models, with limited attention given to improving
the quality of the data. While a few studies have utilized generative
adversarial networks (GANs) to enhance data quality by improving lighting
conditions [22] or removing artifacts from glasses [34], these niche
approaches have demonstrated limited practical success, as they do not address
the fundamental challenge of obtaining high-quality datasets.

###  2.3 Super-Resolution

The first application of GANs to super-resolution was in SRGAN, which
outperformed prior methods and achieved state-of-the-art results [23]. The
authors attributed this, in part, to their use of a perceptual loss function
that accounted for perceptual similarity instead of just similarity in pixel
space. At that time, a common problem with super-resolution was the presence
of artifacts when upsampling. ESRGAN addressed this issue by identifying that
batch normalization layers tended to create unwanted artifacts [45]. They also
improved the perceptual loss function and used Residual-in-Residual Dense
Blocks to generate more realistic images consistently. Additionally, they
later proposed REAL-ESRGAN, which incorporated a u-net discriminator and
spectral normalization [44]. They also developed a complex degradation method
that used only synthetic data for Real-World Degradations.

Blind face restoration is often challenging because it requires prior
knowledge, such as facial geometry, to restore realistic details. GFP-GAN
utilized generative facial prior in the image restoration process and achieved
realistic details and state-of-the-art results [43]. Other attempts at image
restoration have employed transformers. For example, one study utilized a
shifted window transformer as a deep feature extractor in a model called
SwinIR, which achieved state-of-the-art results with up to 67% fewer
parameters [24]. However, GANs for super-resolution are becoming less popular,
as the current state of the art is based on an iterative refinement method
[35], which significantly outperforms previous works.

While most existing methods typically employ classical degradation methods
such as downsampling to generate low-resolution images, Real-ESRGAN utilized
synthetic data and a complex degradation model that aimed to simulate real-
world complex degradations [44]. Their degradation model involved combining
multiple classical degradations, including Gaussian filters for blurring,
downsampling using interpolating methods, adding Gaussian, colour, and other
types of noise, and reducing quality through JPEG compression. However, other
researchers found that the degradation model used by Real-ESRGAN lacked
diversity, and they addressed this issue by expanding the model through random
shuffling of the process, adding different levels of noise and compression,
and introducing processed camera sensor noise and RAW image noise [49].

###  2.4 Self-Supervised Learning

Self-supervised learning has emerged as a powerful technique for learning rich
and meaningful representations from unlabeled data. This technique involves
training a model to predict a pretext task from the input data, which then
results in learning useful representations that can be transferred to
downstream tasks. One popular form of self-supervised learning is contrastive
learning, which aims to learn representations by contrasting positive and
negative pairs of samples [18, 15, 40]. SimCLR [5] has been one of the most
successful and extensively used approaches in this field. It has achieved
state-of-the-art results on various benchmarks, including ImageNet and COCO.
The authors have demonstrated that this technique is effective in pretraining
models for several downstream tasks, such as object detection, instance
segmentation, and semantic segmentation.

##  3 On the Usefulness of SR for Gaze Prediction

###  3.1 Not All SR Models Preserve Gaze

Prior studies have suggested that using face images with higher pixel density
in the eye regions could enhance the accuracy of gaze estimation using deep
learning [2]. Super-resolution refers to the process of increasing the
resolution of an image by recovering or generating high-resolution images from
low-resolution inputs. It has been shown to be useful in various computer
vision tasks, but its usefulness for gaze prediction has not been studied. In
particular, while SR can increase detail and clarity, it is not apparent if
the SR process alters the gaze in the resulting image. Along this line, we
examine two different SR models.

Although several SR models are available, we focused on GFP-GAN [43] and
SwinIR [24] due to their high performance in PSNR and SSIM, which are commonly
used to evaluate the quality of generative image techniques. These models
adopt distinct approaches to SR; GFP-GAN leverages facial priors for image
restoration, while SwinIR utilizes a transformer. Furthermore, GFP-GAN was
pretrained using the FFHQ dataset [20] while SwinIR was pretrained on DIV2K
[19] and Flickr2K [41]. For gaze, we use the Full-Face model [52] which is
widely considered as the basis for much of the current research in appearance-
based gaze tracking. We use the MPIIFaceGaze dataset [51] which comprises
45,000 images of size $448\times 448$, gathered from 15 subjects of different
ethnicities under natural lighting conditions. To evaluate the efficacy of the
proposed method, the images are first downsampled to a size of $112\times 112$
allowing room for 4x upsampling. To replicate real-world degradations, we used
the method used in [49] to degrade the images. Subsequently, we applied one of
the pretrained SR models (GFP-GAN [43] or SwinIR [24]) to upscale the images
back to size $448\times 448$. To ensure that upscaling did not serve as a
confounding variable, we also upscaled the dataset using bicubic interpolation
to create a baseline. We then train and evaluate Full-Face with each of the
three new datasets using a leave-one-out cross-validation approach. This
involved training on 14 of the 15 participants and testing on the one
participant that was excluded. We repeated this process for each participant
to calculate the average point of gaze (POG) error. The epoch with the highest
test performance in Table 1.

The aim of the initial experiment was to assess how various SR techniques
compare to an interpolation baseline when applied to images with and without
degradations. The experiment utilized Full-Face [52] as the gaze model and
downsampled data of size 112.

The results of the scenario without degradations, shown in Table 1 were
unexpected, with GFPGAN performing 8.4% worse than the baseline, while
SwinIR-4x exhibited a 2.1% improvement. In the degraded scenario, all
interpolation methods achieved mediocre results compared to the first
scenario, likely due to the presence of complex degradations. Once again, GFP-
GAN under-performed relative to the baseline, however, SwinIR-4x outperformed
the baseline by 7.3%. The results indicate that while GFP-GAN can generate
visually sharp images, it is inadequate for the task of gaze estimation. On
the other hand, SwinIR outperformed the baseline with both regular and
degraded images. The most interesting finding is that SwinIR achieved a 5.2%
improvement in relative performance indicating that it was more effective in
degraded scenarios, emphasizing the applicability of SR in real-world
situations. While this experiment shows a large contrast in SR methods, it has
also proved that SR can be an effective tool for gaze estimation.

Degradation Type | Interpolation Type | POG
---|---|---
None | Interpolation | 4.20∘
GFP-GAN-4x | 4.59∘
SwinIR-4x | 4.11∘
Complex Degradations | Interpolation | 5.47∘
GFP-GAN-4x | 5.76∘
SwinIR-4x | 5.10∘
Table 1: Comparison of different SR methods for gaze prediction. GFP-GAN
performs worse than a simple interpolation. SwinIR significantly improves the
gaze prediction. In all cases, Full-Face is used as the gaze model.

As demonstrated in Table 1, not all SR methods are suitable for task dependent
image restoration. GFP-GAN’s poor gaze tracking performance might be
attributed to the use of facial priors, which bias the model into
hallucinating facial features from noise. Additionally, ”in the face of
uncertainty”, GFP-GAN’s dependence on facial priors resulted in restored
images having a similar gaze focused on the centre of the screen. This gaze
hallucination is most evident when the input images are significantly
degraded, indicating that GFP-GAN is susceptible to mode collapse. The issue
of mode collapse is apparent in the hallucinated gazes depicted in Figure 2.
Unlike other restoration approaches, GFP-GAN was trained on the FFHQ dataset
[20], which contains mostly images of people looking directly at the camera.
Since GFP-GAN does not attempt to maintain the original gaze and is influenced
by facial priors, the reconstructed faces have their gazes centred and
therefore ”looks can deceiving”.

Figure 2:  Gaze Mode Collapse. Left: original images. Centre: degraded images.
Right: restored images using GFP-GAN. GFP-GAN tends to generate images where
subject starring at the camera. The red arrow depicts the ground truth gaze
vector while the blue is predicted.

###  3.2 Improving Gaze Estimation with SR

In this section, to further improve the gaze prediction results, Full-Face is
replaced with the state-of-the-art model, GazeTR [9] which uses a hybrid CNN
and transformer architecture. Considering the poor performance of GFP-GAN, as
discussed in the previous section, only SwinIR will be used for the remaining
experiments. We demonstrate that the resulting framework surpasses GazeTR and
produces state-of-the-art results on the MPIIFaceGaze dataset [51]. To
reproduce claims by the authors, GazeTR was pretrained on ETH-XGaze which is a
dataset consisting of over one million high-resolution gaze images presented
in multiple head poses. The dataset was collected from 110 subjects using
multiple SLR cameras, lighting conditions and a calibrated ground truth [50].

The proposed method (second section in Table 2) uses additional preprocessing
techniques.“SwinIR-2x Downsampled” upscales the initial $224\times 224$ data
(Table 2 “Input” column) by 2x to a size of $448\times 448$ followed by
downsampling back to $224\times 224$ for a fair comparison to GazeTR. The
$224\times 224$ data is also upsampled to $448\times 448$ using bicubic
interpolation to give a fair comparison. ”SwinIR-2x” uses SR to upscale the
data 2x, and is kept at size $448\times 448$ to compare to the interpolated
baseline mentioned previously. Finally, the original $448\times 448$ high-
resolution data is used with GazeTR as an upper limit on performance.

The first section in Table 2 demonstrates prior works performance on images of
size $224\times 224$, while the second section demonstrates results from our
second experiment. The results show ”SwinIR-2x Downsampled” achieved a POG of
3.94∘ and demonstrated a 1.5% improvement in performance compared to GazeTR’s
4.00∘. The bicubic interpolation $448\times 448$ baseline achieved a POG of
3.99∘ indicating that regardless of input size, GazeTR achieves 4.00∘.
Interestingly, “SwinIR-2x” achieved a POG of 3.90∘ and demonstrated a 2.3%
improvement over the $448\times 448$ baseline, showing greater relative
performance when the super-resolved image is not compressed. Furthermore,
“SwinIR-2x” achieved a 1% improvement in performance over the original
$448\times 448$ high-resolution data (POG of 3.94∘). The results suggest that
SR can demonstrate better performance even with lower-resolution data. These
findings also suggest that denoising and sharpening images play a larger role
than merely increasing the image size. Finally, SR preprocessing achieves
state-of-the-art results and beats GazeTR [9] by 2.3%.

Method | Input | Gaze | POG
---|---|---|---
Full-Face [52] | 224 | 224 | 4.93∘
RT-Gene [12] | 224 | 224 | 4.66∘
Dilated-Net [7] | 224 | 224 | 4.42∘
Gaze360 [21] | 224 | 224 | 4.06∘
Gaze360 [21] | 224 | 224 | 4.06∘
CA-Net [8] | 224 | 224 | 4.27∘
Mnist [51] | 224 | 224 | 6.39∘
GazeNet [53] | 224 | 224 | 5.76∘
GazeTR-Pure [9] | 224 | 224 | 4.74∘
GazeTR-Hybrid [9] | 224 | 224 | 4.00∘
Proposed SwinIR-2x Downsampled | 224 | 224 | 3.94∘
Proposed SwinIR-2x | 224 | 448 | 3.90∘
Table 2: Comparison with prior works. The gaze column denotes the size of the
input into the gaze model.

###  3.3 Low-Resolution and Degraded Images

To further elaborate on the findings from Table 2, the performance of SR on
low-resolution and degraded images was examined. Additionally, while
reproducing results of GazeTR [9], it was discovered that a significant
portion of its performance was attributed to it’s pretraining on the ETH-XGaze
dataset [50]. Thus, the effect of pretraining the gaze model was also
investigated.

To simulate a low-resolution setting, images were downsamped to $112\times
112$ and also $56\times 56$ (8x smaller than the original dataset) allowing
for a comparison between SR and interpolation upscaling. The results of these
experiments are shown in Table 3. To simulate real-world degradations such as
motion blur, pixelation and compression, BSR-GAN’s complex degradation method
[49] was used and again compare SR to interpolation. Furthermore, these three
experiments were repeated for GazeTR without pretraining since pretraining on
a massive dataset is unrealistic and might skew performance. The pretrained
experiments can be found in the second section of Table 3 and Table 4. While
it is not surprising that pretrained models outperformed non-pretrained, SR
exhibited improved performance relative to the interpolation baseline,
irrespective of pretraining, initial image resolution or the addition of
complex degradations. The experiments carried out emphasize the application of
SR as a preprocessing technique, particularly in scenarios where facial images
are low-resolution (Figure 3) or have been affected by real-world degradations
such as motion blur, pixelation, and compression (Figure 1).

Figure 3:  The left is an image downsampled to size 56, the middle is the previous image super-resolved to 224, and the original 224 image. Pretraining | Interpolation Type | Input | Gaze | POG
---|---|---|---|---
| Interpolation | 56 | 224 | 4.81∘
| SwinIR-4x | 56 | 224 | 4.76∘
| Interpolation | 112 | 448 | 4.53∘
No | SwinIR-4x | 112 | 448 | 4.48∘
| Interpolation | 56 | 224 | 4.31∘
| SwinIR-4x | 56 | 224 | 4.22∘
| Interpolation | 112 | 448 | 4.24∘
Yes | SwinIR-4x | 112 | 448 | 4.21∘
Table 3: Results on low-resolution images. GazeTR is used as the gaze model. Pretraining | Interpolation Type | Input | Gaze | POG
---|---|---|---|---
| Interpolation | 112 | 448 | 5.40∘
No | SwinIR-4x | 112 | 448 | 5.33∘
| Interpolation | 112 | 448 | 5.37∘
Yes | SwinIR-4x | 112 | 448 | 5.20∘
Table 4: Results on degraded images. GazeTR is used as the gaze model.

##  4 Training with Less Labelled Data: a Self-supervised Approach based on SR

Research has previously shown that self-supervised learning can be an
effective method for pretraining a backbone network that can be used for
various downstream tasks such as object detection, instance segmentation, and
semantic segmentation [5]. With this in mind, we examine SR for gaze
prediction through the lens of self-supervision. The SR tasks doesn’t require
any explicit labelling of the data. In particular, we propose to obtain a
backbone network trained for SR on a large unlabelled dataset and subsequently
use it with a relatively simple head trained on small labelled gaze dataset.
The result is a sample efficient gaze model that can achieve competitive
results but with much less labelled training data. This in turn opens up
possibility of building gaze models for scenarios where obtaining a large
labelled gaze dataset is challenging, for example gaze tracking in infants
[13] or older adults [4]. Another area that may benefit from this is gaze
tracking in animals [47, 33, 16].

Figure 4: SuperVision architecture.

The intermediate representations extracted during the SR process might capture
useful patterns and features for gaze prediction. Therefore, we propose a new
architecture called “SuperVision” by adding residual connections from earlier
layers in the SR module to later layers in the Gaze model. Concretely, the
SwinIR-4x model with a ResNet18 are concatenated to produce an end-to-end
model for appearance-based gaze tracking. The SR-Gaze model may lose auxiliary
information during the compression of the SR output to an image. This
bottleneck, coupled with the model’s depth, justifies utilizing a residual
connection from earlier layers in the SR module to later layers in the Gaze
model.

Interpolation Type | POG
---|---
Interpolation | 6.26∘ | 6.06∘ | 6.04∘
SwinIR-4x Downsampled | 6.20∘ | 6.01∘ | 5.91∘
SuperVision | 6.17∘ | 5.90∘ | 4.54∘
Data Used | 5% | 10% | 20%
*For reference, POG of GazeTR with 100% of the data is 5.37∘
Table 5:  Training with a small portion of the training data. A simple
ResNet18 is used as the gaze model. The proposed SuperVision method uses 5x
less labeled data and yet outperforms the state-of-the-art method of GazeTR
which uses 100% of training data.

To evaluate the effectiveness of the proposed method, we run experiments and
compare it with SR preprocessing and an interpolation baseline. For this
experiment low-resolution images of size 112 were used so that the SwinIR’s 4x
upsampling outputs a reasonable $448\times 448$ size. These images are then
downsampled as ResNet18 accepts images of size 224. We run experiments with
5%, 10% and 20% of the labelled data. As expected, results in Table 5 show a
trend of improvement as more labelled training data are used. Particularly
noteworthy is the rapid acceleration of relative improvement over the
interpolation baseline, especially for SuperVision, which achieved a POG of
4.54 and demonstrated a 33% improvement at 20% training data. The SuperVision
model was also evaluated without a residual connection, and the results
indicated a 4.6% improvement with the inclusion of a residual connection.
Moreover, we can directly compare SuperVision’s results to those of Table 4
since both models use the same data of 112 resolution with complex
degradations. Despite using only 20% of the data, SuperVision displayed a 15%
improvement over the state-of-the-art model of GazeTR which is pretrained on
the full ETH-XGaze dataset and also uses 100% of the labeled data in
MPIIFaceGaze dataset. These experiments demonstrate SR can be used to achieve
a high performance with only a fraction of labelled data.

##  5 Limitations

One of the primary challenges we encountered during the training phase was
hardware limitations. Due to these constraints, we had to rely on smaller SR
models such as SwinIR-2x instead of SwinIR-4x for several experiments. The
reason being that scaling up from $224\times 224$ by 4x would not have been
feasible on our Tesla V100 GPUs.

##  6 Conclusion

We explored the potential of SR as a preprocessing step for appearance-based
gaze prediction. The study demonstrated that not all SR models are effective
in preserving gaze direction. However, the proposed two-step framework based
on SR achieved state-of-the-art results on the MPIIFaceGaze dataset. The
method was also evaluated on low-resolution and degraded images which
confirmed the effectiveness of SR. Moreover, the newly introduced SuperVision
architecture, which combines SR and ResNet with skip connections, showed
impressive results using only 20% of the labeled data, outperforming the
state-of-the-art GazeTR method which uses 100% of the data by 15%. These
findings demonstrate the potential of the proposed SuperVision method as a
more efficient and effective approach to appearance-based gaze prediction.

## References

  * [1] Kenneth Alberto Funes Mora and Jean-Marc Odobez.  Geometric generative gaze estimation (g3e) for remote rgb-d cameras.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 1773–1780, 2014.
  * [2] Abid Ali and Yong-Guk Kim.  Deep fusion for 3d gaze estimation from natural face images using multi-stream cnns.  IEEE Access, 8:69212–69221, 2020.
  * [3] Yiwei Bao, Yihua Cheng, Yunfei Liu, and Feng Lu.  Adaptive feature fusion network for gaze tracking in mobile tablets.  In 2020 25th International Conference on Pattern Recognition (ICPR), pages 9936–9943. IEEE, 2021.
  * [4] GJ Chapman and Mark Andrew Hollands.  Evidence for a link between changes to gaze behaviour and risk of falling in older adults during adaptive locomotion.  Gait & posture, 24(3):288–294, 2006.
  * [5] Ting Chen, Simon Kornblith, Mohammad Norouzi, and Geoffrey Hinton.  A simple framework for contrastive learning of visual representations.  In Hal Daumé III and Aarti Singh, editors, Proceedings of the 37th International Conference on Machine Learning, volume 119 of Proceedings of Machine Learning Research, pages 1597–1607. PMLR, 13–18 Jul 2020\.
  * [6] Zhaokang Chen and Bertram E Shi.  Appearance-based gaze estimation using dilated-convolutions.  In Asian Conference on Computer Vision, pages 309–324. Springer, 2018.
  * [7] Zhaokang Chen and Bertram E Shi.  Appearance-based gaze estimation using dilated-convolutions.  In Computer Vision–ACCV 2018: 14th Asian Conference on Computer Vision, Perth, Australia, December 2–6, 2018, Revised Selected Papers, Part VI, pages 309–324. Springer, 2019.
  * [8] Yihua Cheng, Shiyao Huang, Fei Wang, Chen Qian, and Feng Lu.  A coarse-to-fine adaptive network for appearance-based gaze estimation.  In Proceedings of the AAAI Conference on Artificial Intelligence, number 07, pages 10623–10630, 2020.
  * [9] Yihua Cheng and Feng Lu.  Gaze estimation using transformer.  arXiv preprint arXiv:2105.14424, 2021.
  * [10] Yihua Cheng, Haofei Wang, Yiwei Bao, and Feng Lu.  Appearance-based gaze estimation with deep learning: A review and benchmark.  arXiv preprint arXiv:2104.12668, 2021.
  * [11] Tobias Fischer, Hyung Jin Chang, and Yiannis Demiris.  Rt-gene: Real-time eye gaze estimation in natural environments.  In Proceedings of the European Conference on Computer Vision (ECCV), September 2018.
  * [12] Tobias Fischer, Hyung Jin Chang, and Yiannis Demiris.  Rt-gene: Real-time eye gaze estimation in natural environments.  In Proceedings of the European conference on computer vision (ECCV), pages 334–352, 2018.
  * [13] John M Franchak, David J Heeger, Uri Hasson, and Karen E Adolph.  Free viewing gaze behavior in infants and adults.  Infancy, 21(3):262–287, 2016.
  * [14] Kenneth Alberto Funes Mora, Florent Monay, and Jean-Marc Odobez.  Eyediap: A database for the development and evaluation of gaze estimation algorithms from rgb and rgb-d cameras.  In Proceedings of the Symposium on Eye Tracking Research and Applications, pages 255–258, 2014.
  * [15] Jean-Bastien Grill, Florian Strub, Florent Altché, Corentin Tallec, Pierre H. Richemond, Elena Buchatskaya, Carl Doersch, Bernardo Ávila Pires, Zhaohan Daniel Guo, Mohammad Gheshlaghi Azar, Bilal Piot, Koray Kavukcuoglu, Rémi Munos, and Michal Valko.  Bootstrap your own latent: A new approach to self-supervised learning.  CoRR, abs/2006.07733, 2020.
  * [16] Kun Guo, Robert G Robertson, Sasan Mahmoodi, Yoav Tadmor, and Malcolm P Young.  How do monkeys view faces?—a study of eye movements.  Experimental Brain Research, 150:363–374, 2003.
  * [17] Dan Witzner Hansen and Qiang Ji.  In the eye of the beholder: A survey of models for eyes and gaze.  IEEE transactions on pattern analysis and machine intelligence, 32(3):478–500, 2009.
  * [18] Kaiming He, Haoqi Fan, Yuxin Wu, Saining Xie, and Ross B. Girshick.  Momentum contrast for unsupervised visual representation learning.  CoRR, abs/1911.05722, 2019.
  * [19] Andrey Ignatov, Radu Timofte, et al.  Pirm challenge on perceptual image enhancement on smartphones: report.  In European Conference on Computer Vision (ECCV) Workshops, January 2019.
  * [20] Tero Karras, Samuli Laine, and Timo Aila.  A style-based generator architecture for generative adversarial networks.  In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, pages 4401–4410, 2019.
  * [21] Petr Kellnhofer, Adria Recasens, Simon Stent, Wojciech Matusik, and Antonio Torralba.  Gaze360: Physically unconstrained gaze estimation in the wild.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 6912–6921, 2019.
  * [22] Jung-Hwa Kim and Jin-Woo Jeong.  Gaze estimation in the dark with generative adversarial networks.  In ACM Symposium on Eye Tracking Research and Applications, pages 1–3, 2020.
  * [23] Christian Ledig, Lucas Theis, Ferenc Huszár, Jose Caballero, Andrew Cunningham, Alejandro Acosta, Andrew Aitken, Alykhan Tejani, Johannes Totz, Zehan Wang, et al.  Photo-realistic single image super-resolution using a generative adversarial network.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 4681–4690, 2017.
  * [24] Jingyun Liang, Jiezhang Cao, Guolei Sun, Kai Zhang, Luc Van Gool, and Radu Timofte.  Swinir: Image restoration using swin transformer.  In Proceedings of the IEEE/CVF international conference on computer vision, pages 1833–1844, 2021.
  * [25] Feng Ling and Giles M Foody.  Super-resolution land cover mapping by deep learning.  Remote Sensing Letters, 10(6):598–606, 2019.
  * [26] Hongying Liu, Zhubo Ruan, Peng Zhao, Chao Dong, Fanhua Shang, Yuanyuan Liu, Linlin Yang, and Radu Timofte.  Video super-resolution based on deep learning: a comprehensive survey.  Artificial Intelligence Review, 55(8):5981–6035, 2022.
  * [27] Feng Lu, Takahiro Okabe, Yusuke Sugano, and Yoichi Sato.  Learning gaze biases with head motion for head pose-free gaze estimation.  Image and Vision Computing, 32(3):169–179, 2014.
  * [28] Feng Lu, Yusuke Sugano, Takahiro Okabe, and Yoichi Sato.  Adaptive linear regression for appearance-based gaze estimation.  IEEE transactions on pattern analysis and machine intelligence, 36(10):2033–2046, 2014.
  * [29] Feng Lu, Yusuke Sugano, Takahiro Okabe, and Yoichi Sato.  Gaze estimation from eye appearance: A head pose-free method via eye image synthesis.  IEEE Transactions on Image Processing, 24(11):3680–3693, 2015.
  * [30] Kenneth Alberto Funes Mora and Jean-Marc Odobez.  Person independent 3d gaze estimation from remote rgb-d cameras.  In 2013 IEEE International Conference on Image Processing, pages 2787–2791. IEEE, 2013.
  * [31] Carlos Hitoshi Morimoto, Arnon Amir, and Myron Flickner.  Detecting eye position and gaze from a single camera and 2 light sources.  In 2002 International Conference on Pattern Recognition, volume 4, pages 314–317. IEEE, 2002.
  * [32] Atsushi Nakazawa and Christian Nitschke.  Point of gaze estimation through corneal surface reflection in an active illumination environment.  In European Conference on Computer Vision, pages 159–172. Springer, 2012.
  * [33] Tadatoshi Ogura, Mizuki Maki, Saki Nagata, and Sanae Nakamura.  Dogs (canis familiaris) gaze at our hands: A preliminary eye-tracker experiment on selective attention in dogs.  Animals, 10(5):755, 2020.
  * [34] Akshay Rangesh, Bowen Zhang, and Mohan M Trivedi.  Driver gaze estimation in the real world: Overcoming the eyeglass challenge.  In 2020 IEEE Intelligent Vehicles Symposium (IV), pages 1054–1059. IEEE, 2020.
  * [35] Chitwan Saharia, Jonathan Ho, William Chan, Tim Salimans, David J. Fleet, and Mohammad Norouzi.  Image super-resolution via iterative refinement, 2021.
  * [36] Wenzhe Shi, Jose Caballero, Christian Ledig, Xiahai Zhuang, Wenjia Bai, Kanwal Bhatia, Antonio M Simoes Monteiro de Marvao, Tim Dawes, Declan O’Regan, and Daniel Rueckert.  Cardiac image super-resolution with global correspondence using multi-atlas patchmatch.  In Medical Image Computing and Computer-Assisted Intervention–MICCAI 2013: 16th International Conference, Nagoya, Japan, September 22-26, 2013, Proceedings, Part III 16, pages 9–16. Springer, 2013\.
  * [37] Yusuke Sugano, Yasuyuki Matsushita, and Yoichi Sato.  Learning-by-synthesis for appearance-based 3d gaze estimation.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 1821–1828, 2014.
  * [38] Yusuke Sugano, Yasuyuki Matsushita, Yoichi Sato, and Hideki Koike.  An incremental learning method for unconstrained gaze estimation.  In European conference on computer vision, pages 656–667. Springer, 2008.
  * [39] Li Sun, Zicheng Liu, and Ming-Ting Sun.  Real time gaze estimation with a consumer depth camera.  Information Sciences, 320:346–360, 2015.
  * [40] Yonglong Tian, Chen Sun, Ben Poole, Dilip Krishnan, Cordelia Schmid, and Phillip Isola.  What makes for good views for contrastive learning.  CoRR, abs/2005.10243, 2020.
  * [41] Radu Timofte, Eirikur Agustsson, Luc Van Gool, Ming-Hsuan Yang, and Lei Zhang.  Ntire 2017 challenge on single image super-resolution: Methods and results.  In Proceedings of the IEEE conference on computer vision and pattern recognition workshops, pages 114–125, 2017.
  * [42] Roberto Valenti, Nicu Sebe, and Theo Gevers.  Combining head pose and eye location information for gaze estimation.  IEEE Transactions on Image Processing, 21(2):802–815, 2011.
  * [43] Xintao Wang, Yu Li, Honglun Zhang, and Ying Shan.  Towards real-world blind face restoration with generative facial prior.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 9168–9178, 2021.
  * [44] Xintao Wang, Liangbin Xie, Chao Dong, and Ying Shan.  Real-esrgan: Training real-world blind super-resolution with pure synthetic data.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 1905–1914, 2021.
  * [45] Xintao Wang, Ke Yu, Shixiang Wu, Jinjin Gu, Yihao Liu, Chao Dong, Yu Qiao, and Chen Change Loy.  Esrgan: Enhanced super-resolution generative adversarial networks.  In Proceedings of the European conference on computer vision (ECCV) workshops, pages 0–0, 2018.
  * [46] Oliver Williams, Andrew Blake, and Roberto Cipolla.  Sparse and semi-supervised visual mapping with the s^ 3gp.  In 2006 IEEE Computer Society Conference on Computer Vision and Pattern Recognition (CVPR’06), volume 1, pages 230–237. IEEE, 2006.
  * [47] Alexander B Wiltschko, Matthew J Johnson, Giuliano Iurilli, Ralph E Peterson, Jesse M Katon, Stan L Pashkovski, Victoria E Abraira, Ryan P Adams, and Sandeep Robert Datta.  Mapping sub-second structure in mouse behavior.  Neuron, 88(6):1121–1135, 2015.
  * [48] Xuehan Xiong, Zicheng Liu, Qin Cai, and Zhengyou Zhang.  Eye gaze tracking using an rgbd camera: A comparison with a rgb solution.  In Proceedings of the 2014 ACM International Joint Conference on Pervasive and Ubiquitous Computing: Adjunct Publication, pages 1113–1121, 2014\.
  * [49] Kai Zhang, Jingyun Liang, Luc Van Gool, and Radu Timofte.  Designing a practical degradation model for deep blind image super-resolution.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 4791–4800, 2021.
  * [50] Xucong Zhang, Seonwook Park, Thabo Beeler, Derek Bradley, Siyu Tang, and Otmar Hilliges.  Eth-xgaze: A large scale dataset for gaze estimation under extreme head pose and gaze variation.  In European Conference on Computer Vision, pages 365–381. Springer, 2020.
  * [51] Xucong Zhang, Yusuke Sugano, Mario Fritz, and Andreas Bulling.  Appearance-based gaze estimation in the wild.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 4511–4520, 2015.
  * [52] Xucong Zhang, Yusuke Sugano, Mario Fritz, and Andreas Bulling.  It’s written all over your face: Full-face appearance-based gaze estimation.  In Proceedings of the IEEE conference on computer vision and pattern recognition workshops, pages 51–60, 2017.
  * [53] Xucong Zhang, Yusuke Sugano, Mario Fritz, and Andreas Bulling.  Mpiigaze: Real-world dataset and deep appearance-based gaze estimation.  IEEE transactions on pattern analysis and machine intelligence, 41(1):162–175, 2017.
  * [54] Zhiwei Zhu, Qiang Ji, and Kristin P Bennett.  Nonlinear eye gaze mapping function estimation via support vector regression.  In 18th International Conference on Pattern Recognition (ICPR’06), volume 1, pages 1132–1135. IEEE, 2006.