Corresponding author: Tanay Bhandarkar
E-mail<EMAIL_ADDRESS>
# The Simons Observatory: Development and Validation of the Large Aperture
Telescope Receiver

Tanay Bhandarkar  Department of Physics and Astronomy, University of
Pennsylvania, Philadelphia, PA, USA  Sanah Bhimani  Department of Physics,
Yale University, New Haven, CT, USA  Gabriele Coppi  Department of Physics,
University of Milano-Bicocca, Milano, Italy  National Institute for Nuclear
Physics , Sezione di Milano-Bicocca, Milano, Italy  Simon Dicker  Department
of Physics and Astronomy, University of Pennsylvania, Philadelphia, PA, USA
Saianeesh K. Haridas  Department of Physics and Astronomy, University of
Pennsylvania, Philadelphia, PA, USA  Kathleen Harrington  Department of
Astronomy and Astrophysics, University of Chicago, Chicago, IL, USA  Jeffrey
Iuliano  Department of Physics and Astronomy, University of Pennsylvania,
Philadelphia, PA, USA  Bradley Johnson  Department of Astronomy, University of
Virginia, Charlottesville, VA, USA  Anna M. Kofman  Department of Physics and
Astronomy, University of Pennsylvania, Philadelphia, PA, USA  Jack Lashner
Department of Physics and Astronomy, University of Southern California, LA,
CA, USA  Jenna Moore  School of Earth and Space Exploration, Arizona State
University, Tempe, AZ, USA  David V. Nguyen  Department of Physics, Yale
University, New Haven, CT, USA  John Orlowski-Scherer  Department of Physics
and Astronomy, University of Pennsylvania, Philadelphia, PA, USA  Karen Perez
Sarmiento  Department of Physics and Astronomy, University of Pennsylvania,
Philadelphia, PA, USA  Julia Robe  Department of Physics and Astronomy,
University of Pennsylvania, Philadelphia, PA, USA  Maximiliano Silva-Feaver
Department of Physics, University of California San Diego, La Jolla, CA, USA
Robert J. Thornton  Department of Physics, West Chester University of
Pennsylvania, West Chester, PA, USA  Department of Physics and Astronomy,
University of Pennsylvania, Philadelphia, PA, USA  Yuhan Wang  Department of
Physics, Princeton University, Princeton, NJ, USA  Zhilei Xu  MIT Kavli
Institute, Massachusetts Institute of Technology, Cambridge, MA, USA

###### Abstract

The Simons Observatory (SO) is a ground-based cosmic microwave background
(CMB) survey experiment that consists of three 0.5 m small-aperture telescopes
(SATs) and one 6 m large-aperture telescope (LAT), sited at an elevation of
5200 m in the Atacama Desert in Chile. In order to meet the sensitivity
requirements set for next-generation CMB telescopes, the LAT will deploy
30,000 transition edge sensor (TES) detectors at 100 mK across 7 optics tubes
(OT), all within the Large Aperture Telescope Receiver (LATR). Additionally,
the LATR has the capability to expand to 62,000 TES across 13 OTs. The LAT
will be capable of making arcminute-resolution observations of the CMB, with
detector bands centered at 30, 40, 90, 150, 230, and 280 GHz. We have
rigorously tested the LATR systems prior to deployment in order to fully
characterize the instrument and show that it can achieve the desired
sensitivity levels. We show that the LATR meets cryogenic and mechanical
requirements, and maintains acceptably low baseline readout noise.

###### keywords:

Astronomical Instrumentation, Cosmic Microwave Background, Cryogenic Receiver,
Multiplexing readout, Observational Cosmology

##  1 Introduction

Observations of the cosmic microwave background (CMB) over the last several
decades have been instrumental to the development of modern cosmological
models, such as $\Lambda$-CDM. Our understanding of the CMB has evolved
through full-sky surveys conducted by space-based satellite experiments such
as COBE, WMAP[1] and Planck[2]. Large aperture ground based experiments, such
as the Atacama Cosmology Telescope (ACT)[3] and the South Pole Telescope
(SPT)[4], have improved our ability to make small angular measurements of the
CMB. These small angular measurements have offered insights into the
development of large scale structure, the upper bound on the sum of neutrino
masses, the Sunyaev-Zeldovich effects (SZ) and primordial power spectra. Other
ground based telescopes, such as CLASS[5] and BICEP/Keck, utilize small
apertures to primarily measure the B-mode polarization signal.

The Simons Observatory (SO) is a next-generation, ground-based CMB experiment,
consisting of three 0.5 m small-aperture telescopes (SATs) and one 6 m large-
aperture telescope (LAT). These telescopes will be deployed at an elevation of
5200 m in the Atacama Desert of Chile. The LAT will utilize a crossed-Dragone
optical design, enabling a large focal plane for the telescope. The Large
Aperture Telescope Receiver (LATR) will occupy the LAT’s focal plane with
$\sim$30,000 transition edge sensor (TES) bolometers across 7 optics tubes
(OT), and has the capability to upgrade to $\sim$62,000 TES across 13 OTs.
Each OT is specifically designed for a separate frequency band (with targeted
band centers); low-frequency (30 and 40 GHz), medium-frequency (90 and 150
GHz) and ultra-high frequency (230 and 280 GHz). Over the course of it’s
observing campaign, the LAT will scan 40% of the sky at arcminute resolution,
with significant overlap with other surveys at complementary wavelengths, such
as the Dark Energy Survey (DES)[6], the Dark Energy Spectroscopic Instrument
(DESI)[7] and the Legacy Survey of Space and Time (LSST) with the Vera Rubin
Observatory[8] (Figure 1). The field of view of the LATR will be approximately
$7.05^{\circ}$ in diameter. The LAT will focus on small-angular sciences,
including measurements of the temperature, polarization and lensing power
spectra, the detection of galaxy clusters through the thermal SZ effect, the
detection of extra-galactic point sources and transient event monitoring [9].

Figure 1: Planned sky coverage for the LAT [9].

##  2 LATR Design

The LATR consists of 5 separate nominal temperature stages; 300 K, 80 K, 40 K,
4 K, 1 K and 100 mK. We will describe the design of the cryostat in this
section, and further details can be found in other publications.[10, 11].
Figure 3 shows the cross-sectional view of a full LATR.

###  2.1 Cryo-mechanical and Optical Design

The 300 K vacuum shell is 2.4 m in diameter, and 2.6 m in length, made of
aluminum-6061. The front plate is 6 cm thick with a honeycomb pattern for 13
anti-reflective (AR) coated ultra-high-molecular-weight polyethylene (UHMWPE)
1/8” thick windows. The thicknesses of the vacuum shell, front and back
plates, and windows were determined following a finite element analysis (FEA)
study into the structural requirements at atmospheric pressures [12]. Behind
each UHMWPE window sits a double-sided IR blocker (DSIR), which reflects
thermal radiation away from the the colder stages [13]. The 300 K stage will
mate the LATR to the LAT through a pair of co-rotating rails, that turn the
receiver as the elevation structure of the telescope rotates (see Figure 2).

Figure 2: Computer aided design models of the LAT (left) and LATR (right). The
LATR will sit on the co-rotator rails inside the LAT’s receiver cabin,
enabling it to rotate along with the LAT’s elevation structure. The LAT sits
on a large conic base that will rotate the telescope azimuthally. Figure from
Xu et al. 2021.[11]

The 80 K stage consists of the 80 K filter plate and the 80 K shell. This
filter plate holds a DSIR and an anti-reflection (AR) coated alumina
filter[14] per OT. This set of filters absorb a significant fraction of
radiative power that comes in through the 300 K windows. The 40 K stage holds
a single 40 K DSIR per OT. The 4 K stage consists of a back lid, and filter
plate and a radiation shield. Because the OTs mount on this stage, the 4 K
filter plate is significantly thicker than the other cold filter plates to
better support the weight of the OTs and improve the thermal path from the
pulse tubes to the OTs. Additionally, the 80 K, 40 K and 4 K stages are all
wrapped in 30 layers of multi-layer insulation (MLI), preventing thermal
radiation from reaching the colder stages.

The 80 K stage is cooled by two PT-90 pulse tube coolers, while the 40 K and 4
K stages are cooled by a pair of two-stage PT-420 pulse tube
coolers111https://www.cryomech.com/. All pulse tube heads are coupled to their
respective cold stage by a set of copper straps. Each copper strap is a pair
of oxygen-free, high conductivity (OFHC) copper plates, bridged by copper
braids.

The requirements for the LATR design demands the ability to balance mechanical
integrity with thermal isolation between the different temperature stages. To
accomplish this, the LATR utilizes three separate rings of G-10 tabs. These
tabs mechanically connect the 300 K stage to the 80 K stage, the 300 K stage
to the 40 K stage, and the 40 K stage to the 4 K stage. Each G-10 tab consists
of two aluminum ”feet” epoxyed to a sheet of G-10. The feet are then bolted to
their respective cold stage cold stages.

The TES detectors are cooled to 100 mK in order to reduce thermal noise and
improve sensitivity. To avoid larger thermal gradients across OT 100 mK
stages, the LATR implements a wheel-like structure for its 1 K and 100 mK
thermal back-up structure (BUS) cold stage (see Figure 4). The thermal BUS
stages are both made of OFHC copper. Thermal FEA predicts a gradient of $\sim$
5 mK across the entire 100 mK BUS. These cold stages are cooled with a
Bluefors222https://bluefors.com/ model LD 400 dilution refrigerator (DR). The
thermal BUS stages are coupled to their respective stages of the DR through
gold-plated, braided copper straps. In order to cool the 1 K stage and focal
plane baseplates (FPB) that supports the detector arrays in each OT, the LATR
utilizes a pair copper rod cold fingers (one on the FPB(100 mK) or radiation
lid (1 K), one on each stage of the thermal BUS), bridged by a separate gold-
plated, braided copper strap. The 1 K BUS is mounted to the 4 K plate by six
twill-ply carbon fiber tube tripods. The 100 mK BUS is supported by a truss of
carbon fiber connecting to the 1 K BUS.

Figure 3: A cross-sectional view of the LATR, with the 300 K (gray), 80 K
(gray), 40 K (green), 4 K (pink) stages all visible. Additionally, the 1 K BUS
(brown) and 100 mK BUS (indigo) can be seen at the back of the LATR, close to
the DR. This view also displays the entire filter stack in an optics tube. On
the exterior of the LATR, the SMuRF readout crates can be seen mounted at
different locations, as well as the DR’s turbo pump station on top. [10].

###  2.2 Readout and Detector Design

In order to achieve the sensitivity goals of SO, the LATR will field
$\sim$62,000 TES detectors in its final configuration. To efficiently read out
this large volume of sensors, the LATR implements a microwave squid
multiplexing readout system, carrying tones along coaxial lines from the 300 K
stage down to the the 100 mK universal focalplane modules (UFMs), and back
out. The expectation is to be able to read out $\mathcal{O}(10^{3})$ detectors
per coaxial line.[15] UFM design has been detailed in prior publications.[16]

####  2.2.1 Cold Readout design

The cold readout chain carries RF (4-6 GHz) tones and DC ($<$1 MHz) biases
output by warm readout electronics, into the cryogenic receiver, to the
cryogenic detectors and back out to the SMuRF. From 300 K to 4 K, the LATR
readout system is contained within the Universal Readout Harness (URH), an
assembly that is used throughout SO receivers[17]. In a single URH, there are
up to 12 bias channels and 24 coaxial channel. A bias channel is a 50 pin
input that carries 2 flux ramp line pairs, 12 detector bias pair lines, and 4
cryogenic amplifier bias lines. Each coaxial channel consists of an input
connection and output connection, and the respective cold readout components.
Ultimately, each UFM is read out through two coaxial chains, and each bias
channel in the URH is coupled to two separate coaxial chains. In sum, a
single, fully populated LATR URH can read out up to 12 UFMs.

At the 4 K stage, a series of coaxial connections (referred to as coax
highways from here out) and DC bias lines connect the URHs to individual OTs.
The input lines in the OTs consist of a series of attenuators to optimize tone
power that reaches the resonators, and DC blocks that provide thermal
isolation in coaxial connections between cold stages. On the output line from
the UFMS, there are low-noise amplifiers (LNA) to amplify the output signal,
and isolators to prevent standing waves between the UFM and the LNAs. In
parallel to the coaxial connections, each readout chain contains DC bias lines
that travel all the way down to the UFM. These bias lines provide voltage
biases for the LNAs, detector biases and a ramp signal for the flux ramp
signal used by warm electronics to linearize the resonator frequency shifts
produced by the detectors.[18]

####  2.2.2 Warm Readout design

Each TES on an UFM is coupled to a resonator at a unique resonant frequency,
between 4-6 GHz. SO will uses the SLAC Microresonator Radio Frequency (SMuRF)
electronics to generate these 4-6 GHz frequency tones, transmitted through the
cryogenic RF components via coaxial lines. Up to 147 TES detectors share the
same voltage bias line and 882 detector channels share the same resonator flux
ramping lines, both of which are generated by the SMuRF electronics. [19]

4-6 GHz RF tones are generated in two Advanced Mezzanine Cards (AMC) seated in
a carrier card. Every AMC is connected to a single URH coaxial channel, along
with a warm amplifier mounted on the URH. Each carrier card is coupled to a
Rear Transition Module (RTM), which is responsible for generating all the low
frequency signals (detector biases, flux ramps, amplifier biases). The RTMs
are connected to a single cryostat card (cryocard), which is also mounted
above the URH. The cryocards convert the TES bias and flux ramp signals from
the RTM’s low impedance voltage sources to high impedance current sources. It
also provides regulated power supplies, bias line filtering and gate bias
voltages to LNAs.

In order to reduce the effects of ambient RF noise and minimize cable loss,
the SMuRF crates are mounted directly on the LATR shell. There are extensive
ongoing studies into effects of thermal and vibrational environments on the
warm readout system.

##  3 Validation

###  3.1 Cryogenic Validation

In order to meet detector sensitivity and noise requirements, it is crucial
for the LATR to cool every OT FPB to a stable 100 mK level. Our 100 mK loading
estimates for various configurations are derived from baseline load curves. In
other words, we first characterized our DR by applying a range of thermal
loads to its mixing chamber and still line stages. We then methodically built
up the LATR by adding different cryogenic components; the thermal BUS, dark
OTs, and finally optically coupled OTs. Past studies of the LATR have
calculated the thermal load of the 100 mK thermal BUS to be $\sim$20$\mu$W,
and the thermal load of a single dark OT to be$\lesssim$5$\mu$W, and the
temperature of the OT’s FPB varied by less than 0.3 mK over 24 hours [10] (In
this study, ”dark” indicates no cold optics below 40 K, and an aluminum blank
over the OT’s 4 K filter cells). Here we will detail the LATR cryogenic
thermal loading values with 7 dark OTs, and then with 1 optical. Details on 80
K, 40 K and 4 K loading can be found in Tables 1,2 and 3.

####  3.1.1 7 Dark OTs

Figure 4: The LATR with its back plates removed. seven optics tubes (OTs) are
installed in this photo. Each OT is connected to one of the two universal
readout harnesses (URH) via a lengthy coaxial connection, colloquially named
the ”Coax highway” at 4 K.

A significant milestone in LATR lab testing was the cryo-mechanical
integration of 7 dark OTs, in the deployment configuration positions, as shown
in Figure 5. In this configuration, the OTs were covered with an aluminum
plate at the 4 K filter stack. At the 100 mK stage, the DR mixing chamber
plate stabilized at about 43.9 mK. Based on prior load characterizations, we
estimated that each OT on average exerts $\sim$7.1$\mu$W at 100 mK through
parasitic loading. The hottest FPB in this configuration was $\sim$75 mK.
Projecting these numbers out to a 13 OT configuration, we estimate that the
hottest FPB would be $\leq$ 105 mK. At this temperature, we run the risk of
putting detectors over their expected baseline temperature. In order to
rectify this, we plan to install additional thermal straps coupling the DR to
the 100 mK thermal BUS, thus improving the transfer of cooling power to the
BUS.

At the 40 K and 4 K stage, we estimate $50.9\pm 0.9$W and $1.27\pm 0.06$W
respectively. Extrapolating this loading to 13 OTs, we find that the 40 K
stage and 4 K stage will cool to $\sim 32.5$K and $3.95$ K respectively. Both
of these values fall within specifications. At 80 K, we measure $87.9$W. For
13 OTs, these loading values suggest at filter plate temperature of $\sim$82
K, which we expect will be sufficiently cold.

When the LATR is finally in the field, cooldown speed is of great importance;
the more time it takes to reach base temperature, the less time is spent
observing. With 7 OTs, we found that it takes approximately 8 days for all of
the LATR cold stages to reach base temperature (Figure 6 ). The cooldown time
is dominated by the drop to 4 K on the coldest stages, while the 80 K and 40 K
stages reach their respective base temperatures within 5 days. We utilize a
mechanical heat switch assembly, which creates thermal shorts between the 4 K
plate and the 1 K and 100 mK BUS. Without this assembly, simulations suggest
it would the LATR to take over 35 days to reach base temperature. [20]. We
expect that with 13 OTs, the heat switches will cut that time in half.

Figure 5: The deployment configuration of the LATR with 7 optics tubes (OTs),
as viewed from the rear of the receiver. The color of the OT indicates which
universal readout harness (URH) it will be read out through. OTs i1, i3, i4
and i6 will be MF tubes. OT c1 and i5 will be UHF tubes. OTo6 will be the only
LF tube intially deployed. Figure 6: Cooldown times for the LATR in various
configurations. In all of these configurations, the optics tubes (OT) are
capped at 4 K.

####  3.1.2 1 Optical OT

After completing dark cryogenic testing of the OTs, we installed one OT
equipped with the full cold filter and lens stack at the OTi6 position (Figure
5). This test configuration included the full complement of readout hardware,
but did not include full detector array assemblies. In lieu of detector
arrays, we installed three copper mock arrays with Al feedhorns, to mimic the
expected mass and collecting area of an UFM at that stage. The base
temperature of the OT’s FPB was $\sim$55 mK, with very little variation over
the course of 24 hours (see Figure 7). Based on thermal load curves taken when
the LATR was dark, and only cooling the 100 mK thermal BUS, we determined that
the thermal load of a single optical OT $\leq$6.0 $\mu$W. Thermal simulations
predict that the full detector arrays will add $\sim$1.0 $\mu$W per OT when
they are voltage biased. This loading number is higher than predicted at the
100 mK stage; we expect additional straps mentioned in section 3.1.1 to
augment cooling transfer sufficiently.

Based on our cryogenic testing, there appear to be small discrepancies between
loading numbers extracted from different cooldowns. There are a couple
explanations for this phenomenon. One, comparing thermal loading between
cooldowns assumes that cold strapping is assembled identically between tests.
In reality, this is difficult to ensure; the strapping surfaces experience
wear and bolts may not always be tightened to the same torque. Second, it is
difficult to precisely measure the thermal loading to sub-$\mu$W levels due to
variations in hardware and noise in thermometry data.

Figure 7: Temperature of the OT focal plane baseplate (FPB) behind the full optical stack, over the course of 24 hours. During this time, the cryostat was left untouched. The inset plot displays the power spectral density (PSD) of that thermometer during the same period of time. The flatness of the PSD indicates the stability of the temperature. Test Configuration | 80 K Plate Temperature | Thermal Loading
---|---|---
Dark |  $35.4-38.8$ K |  $22_{-1}^{+1}$ W
2 Filter Stacks |  $37.2-43.1$ K |  $35_{-1}^{+1}$ W
3 Filter Stacks |  $41.0-48.7$ K |  $42_{-1}^{+1}$ W
7 Filter Stacks |  $46.8-58.3$ K |  $69_{-1}^{+1}$ W
Table 1: Loading estimates on the 80 K stage, based on temperatures and load curves taken from the two PT-90 cryocoolers. The 80 K Plate temperatures represent the range of temperatures between the warmest point on the filter point to the coldest pulse tube. Test Configuration | 40 K Plate Temperature | 40 K Thermal Loading
---|---|---
Dark* |  $29.3-26.8$ K |  $33_{-1}^{+1}$ W
2 Filter Stacks |  $27.2-48.5$ K |  $34_{-1}^{+1}$ W
3 Filter Stacks |  $27.3-49.1$ K |  $36_{-1}^{+1}$ W
7 Dark OTs |  $30.7-53.7$ K |  $51_{-1}^{+1}$ W
1 Optical OT |  $27.6-49.1$ K |  $37_{-1}^{+1}$ W
Table 2: Loading estimates on the 40 K, based on temperatures and load curves taken from the first stages of the two PT-420 cryocoolers. The asterisk indicates a cooldown with no URH installed. Test configurations with OTs assume at least 7 full filter stacks down through 40 K. Plate temperatures represent the range of temperatures between the warmest point on the filter point to the coldest pulse tube.  Test Configuration | 4 K Plate Temperature | 4 K Thermal Loading
---|---|---
Dark* |  $2.74-5.02$ K |  $0.8_{-0.1}^{+0.1}$ W
2 Filter Stacks |  $3.04-6.68$ K |  $0.8_{-0.1}^{+0.2}$ W
3 Filter Stacks |  $2.90-5.15$ K |  $1.3_{-0.2}^{+0.2}$ W
7 Dark OTs |  $3.09-6.44$ K |  $1.27_{-0.1}^{+0.1}$ W
1 Optical OT |  $2.83-5.15$ K |  $0.8_{-0.1}^{+0.1}$ W
Table 3: Loading estimates on the 4 K stages, based on temperatures and load
curves taken from the second stage of two PT-420 cryocoolers. The asterisk
indicates a cooldown with no URH installed. Test configurations with OTs
assume at least 7 full filter stacks down through 40 K. Plate temperatures
represent the range of temperatures between the warmest point on the filter
point to the coldest pulse tube.

###  3.2 Readout Validation

Given the complexity of the LATR’s readout design, it is important to show
that we are able to maintain acceptably low readout noise when deploying UFMs.
Past studies have shown TES white noise levels are close to expected models in
a single pixel box, containing 6 TES and 65 resonators[11]. LATR MF UFM in
situ testing is planned before deployment. This testing includes establishing
baseline receiver noise levels, and studies into the effects of various
environmental factors, such as vibrations, cryogenic thermal instability,
ambient RF backgrounds and warm readout thermal coupling. Optical testing,
including holography, cross-polarization measurements and beam mapping, has
been conducted in a testbed cryostat [21].

##  4 Conclusion

The Simons Observatory (SO) Large Aperture Telescope (LAT) will be a sensitive
millimeter-wave telescope with a large focal plane. The Large Aperture
Telescope Receiver (LATR) will occupy this focal plane with up to $\sim$62,000
transition edge sensor (TES) detectors, across 13 modular optics tubes (OT).
Each OT will contain 3 universal focal plane module (UFM) arrays, and the
respective readout and cold optics colder than 4 K. Due to the sensitivity
requirements, the LATR must balance a complex optical and readout design with
a robust cryo-mechanical design. Through several rounds of cryogenic testing,
we show that the LATR will be capable of cooling down cold optical elements,
readout assemblies and detector arrays for 7 OTs for its initial deployment.
For 13 OTs, we will add additional cold strapping to ensure we are comfortably
within cryogenic specifications. In order to optimize TES detector noise
levels, the LATR must achieve a stable 100 mK base temperature. Cryogenic lab
testing of the receiver shows that it is capable of cooling the focal plane
baseplate of a single optically coupled optics tube (OT) to a base temperature
of $\sim$54 mK. Testing of multiple optically coupled OTs with fully
functional UFMs is underway.

Additionally, lab testing has demonstrated the viability of the LATR’s readout
system, including its ability to track resonators and bias detectors. Further
laboratory readout testing is planned, including studies into baseline noise
levels and environmental effects. These tests will employ deployment-grade MF
UFMs, with full cold optical stacks installed. At the time of this proceeding,
the expectation is that the LATR will deployed to the field in early 2023,
with first science observations occurring later that year.

##  5 Acknowledgements

This work was funded by the Simons Foundation (Award #457687, B.K.) and the
University of Pennsylvania. Gabriele Coppi is supported by the European
Research Council under the Marie Skłodowska Curie actions through the
Individual European Fellowship No. 892174 PROTOCALC.

## References

  * [1] Bennett, C. L., Larson, D., Weiland, J. L., Jarosik, N., Hinshaw, G., Odegard, N., Smith, K. M., Hill, R. S., Gold, B., Halpern, M., Komatsu, E., Nolta, M. R., Page, L., Spergel, D. N., Wollack, E., Dunkley, J., Kogut, A., Limon, M., Meyer, S. S., Tucker, G. S., and Wright, E. L., “NINE-YEAR WILKINSON MICROWAVE ANISOTROPY PROBE ( WMAP ) OBSERVATIONS: FINAL MAPS AND RESULTS,” The Astrophysical Journal Supplement Series 208, 20 (sep 2013).
  * [2] Planck Collaboration, “Planck 2018 results - i. overview and the cosmological legacy of planck,” A&A 641, A1 (2020).
  * [3] Thornton, R. J. et al., “The atacama cosmology telescope: The polarization-sensitive actpol instrument,” The Astrophysical Journal Supplement Series 227, 21 (Dec 2016).
  * [4] Ruhl, J. et al., “The south pole telescope,” Millimeter and Submillimeter Detectors for Astronomy II (Oct 2004).
  * [5] Xu, Z., Brewer, M. K., Rojas, P. F., Li, Y., Osumi, K., Pradenas, B., et al., “Two-year Cosmology Large Angular Scale Surveyor (CLASS) Observations: 40 GHz Telescope Pointing, Beam Profile, Window Function, and Polarization Performance,” The Astrophysical Journal 891, 134 (Mar 2020).
  * [6] Sevilla-Noarbe, I., Bechtol, K., Kind, M. C., Rosell, A. C., Becker, M. R., Drlica-Wagner, A., Gruendl, R. A., Rykoff, E. S., Sheldon, E., Yanny, B., et al., “Dark energy survey year 3 results: Photometric data set for cosmology,” The Astrophysical Journal Supplement Series 254, 24 (may 2021).
  * [7] DESI Collaboration, “The desi experiment part i: Science,targeting, and survey design,” (2016).
  * [8] Ivezić, Ž., Kahn, S. M., Tyson, J. A., et al., “LSST: From science drivers to reference design and anticipated data products,” The Astrophysical Journal 873, 111 (mar 2019).
  * [9] The Simons Observation Collaboration, “The simons observatory: science goals and forecasts,” Journal of Cosmology and Astroparticle Physics 2019, 056–056 (feb 2019).
  * [10] Zhu, N., Bhandarkar, T., Coppi, G., Kofman, A. M., Orlowski-Scherer, J. L., Xu, Z., et al., “The Simons Observatory Large Aperture Telescope Receiver,” The Astrophysical Journal Supplement 256, 23 (Sept. 2021).
  * [11] Xu, Z. et al., “The simons observatory: The large aperture telescope (LAT),” Research Notes of the AAS 5, 100 (Apr 2021).
  * [12] Orlowski-Scherer, J. L., Zhu, N., Xu, Z., et al., “Simons Observatory large aperture receiver simulation overview,” in [Millimeter, Submillimeter, and Far-Infrared Detectors and Instrumentation for Astronomy IX ], Zmuidzinas, J. and Gao, J.-R., eds., 10708, 644 – 657, International Society for Optics and Photonics, SPIE (2018).
  * [13] Ade, P. A. R., Pisano, G., Tucker, C., and Weaver, S., “A review of metal mesh filters,” in [Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series ], Zmuidzinas, J., Holland, W. S., Withington, S., and Duncan, W. D., eds., Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series 6275, 62750U (June 2006).
  * [14] Golec, J. E., McMahon, J. J., Ali, A., Dicker, S., Galitzki, N., Harrington, K., Westbrook, B., Wollack, E. J., Xu, Z., and Zhu, N., “Design and fabrication of metamaterial anti-reflection coatings for the Simons Observatory,” in [2021SPIE11451 ], Navarro, R. and Geyl, R., eds., 11451, 1182 – 1189, International Society for Optics and Photonics, SPIE (2020).
  * [15] Rao, M. S., Silva-Feaver, M., et al., “Simons Observatory Microwave SQUID Multiplexing Readout – Cryogenic RF Amplifier and Coaxial Chain Design,” arXiv e-prints , arXiv:2003.08949 (Mar. 2020).
  * [16] McCarrick, H. and others, E. H., “The simons observatory microwave SQUID multiplexing detector module design,” The Astrophysical Journal 922, 38 (nov 2021).
  * [17] Moore, J. et al., “Development and performance of Universal Readout Harnesses for the Simons Observatory,” in [Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series ], (2022).  in preparation.
  * [18] Mates, J. A. B., Irwin, K. D., Vale, L. R., Hilton, G. C., Gao, J., and Lehnert, K. W., “Flux-ramp modulation for squid multiplexing,” JLTP 167, 707–712 (Jun 2012).
  * [19] Henderson, S. W. et al., “Highly-multiplexed microwave SQUID readout using the SLAC Microresonator Radio Frequency (SMuRF) electronics for future CMB and sub-millimeter surveys,” in [Millimeter, Submillimeter, and Far-Infrared Detectors and Instrumentation for Astronomy IX ], Zmuidzinas, J. and Gao, J.-R., eds., 10708, 170 – 185, International Society for Optics and Photonics, SPIE (2018).
  * [20] Coppi, G., Xu, Z., et al., “Cooldown strategies and transient thermal simulations for the Simons Observatory,” in [Millimeter, Submillimeter, and Far-Infrared Detectors and Instrumentation for Astronomy IX ], Zmuidzinas, J. and Gao, J.-R., eds., 10708, 246 – 258, International Society for Optics and Photonics, SPIE (2018).
  * [21] Harrington, K. et al., “The Integration and Testing Program for the Simons Observatory Large Aperture Telescope Optics Tubes,” in [Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series ], (Dec. 2020).


.pdf

# FLECS: A Federated Learning Second-Order Framework via Compression and
Sketching

Artem Agafonov  Moscow Institute of Physics and Technology, Dolgoprudny,
Russia  Mohamed bin Zayed University of Artificial Intelligence, Abu Dhabi,
UAE  Dmitry Kamzolov  Mohamed bin Zayed University of Artificial Intelligence,
Abu Dhabi, UAE  Moscow Institute of Physics and Technology, Dolgoprudny,
Russia  Rachael Tappenden  University of Canterbury, Christchurch 8041, New
Zealand  Alexander Gasnikov  Moscow Institute of Physics and Technology,
Dolgoprudny, Russia  Institute for Information Transmission Problems, Moscow,
Russia  HSE University, Moscow, Russia  Martin Takáč  Mohamed bin Zayed
University of Artificial Intelligence, Abu Dhabi, UAE

###### Abstract

Inspired by the recent work FedNL (Safaryan et al, FedNL: Making Newton-Type
Methods Applicable to Federated Learning), we propose a new communication
efficient second-order framework for Federated learning, namely FLECS. The
proposed method reduces the high-memory requirements of FedNL by the usage of
an L-SR1 type update for the Hessian approximation which is stored on the
central server. A low dimensional ‘sketch’ of the Hessian is all that is
needed by each device to generate an update, so that memory costs as well as
number of Hessian-vector products for the agent are low. Biased and unbiased
compressions are utilized to make communication costs also low. Convergence
guarantees for FLECS are provided in both the strongly convex, and nonconvex
cases, and local linear convergence is also established under strong
convexity. Numerical experiments confirm the practical benefits of this new
FLECS algorithm.

##  1 Introduction

This paper presents a novel second-order method for Federated Learning called
FLECS. In this setup, multiple agents/devices/workers collaborate to solve a
machine learning problem by communicating over a central node (server). It is
assumed that the central node is a large and powerful machine, but the
agents/devices are much smaller (for example, a cellphone, or laptop).
Moreover, raw data is stored locally on each worker and cannot be exchanged or
transferred. The finite sum minimization, or Federated Learning, problem can
be formalized in the following way:

$\textstyle{\min}_{w\in\mathbb{R}^{d}}\left\\{F(w):=\tfrac{1}{n}\textstyle{\sum}_{i=1}^{n}f_{i}(w)\right\\},$
(1)

where $f_{i}(w)$ is the loss function (parametrized by $w\in\mathbb{R}^{d}$)
associated with the data stored on the $i$-th machine, $F(w)$ is the empirical
loss function, and both the problem dimension $d$, and the number of machines
$n$, can be large.
Problems of the form (1) arise not only in distributed machine learning [23,
24], but also in robotics, resource allocation, power system control, control
of drone or satellite networks, distributed statistical inference and optimal
transport, and multi-agent reinforcement learning [39, 33, 28, 14, 37, 35].
A variety of first-ordered methods have been introduced for the problem (1)
under the stated computational set-up (a central server connected to multiple
smaller devices). Early work includes that of [10, 11, 26], which focused on
distributed first order optimization algorithms that were designed with the
computational structure in mind. Recent advances include works focusing on
reducing the communication burden via compression techniques [2, 19, 9],
including momentum for algorithm speed up [27, 43], variance reduction
techniques [20], the use of adaptive learning rates [40], and works that focus
on maintaining client data privacy [16].
Beyond first-order methods, second-order methods have also been developed for
Federated Learning. Some of these algorithms [42, 15, 12, 5, 1] utilize
statistical similarity (the homogeneous data setting), which, roughly
speaking, means that the local functions $f_{i}$ are similar, and they
approximate the overall empirical risk function $F$ well. FedNL [36] for the
second-order method for Federated Learning in the case of _heterogeneous_
data. While FedNL makes a valuable step toward the applicability of second-
order methods to Federated Learning, this framework seems to be impractical
for large-scaled problems, since $d\times d$ Hessian approximations are stored
locally on workers and transferred to the server via compressed matrix
communications. (This violates the assumption in the current work, that the
individual workers have small memory.)
The goal here is to build upon the work in [36] and develop a new second order
framework for Federated Learning that is based upon more realistic assumptions
about the available computational architecture. That is, we assume that
individual workers are devices with relatively low memory availability and
computational capacity (e.g. mobile phones), but they can connect to a large
and powerful central machine. The typical situation, that downloading is much
faster than uploading, is also assumed. 111 An example of such a setting is
(e.g., battery powered) wireless network [34, 41] where the cost of sending
data from devices to the centralized server is energy demanding, implying
severely limited communication bandwidth to the server [18, 8]. On the other
hand, the server is usually connected to an energy network and has no
constraints for sending data to other devices in the network. We introduce the
FLECS second-order framework for strongly convex and non-convex optimization,
which places no additional memory on workers and reduces the number Hessian-
vector products per node. Therefore, our algorithm can be viewed as the
lightweight FedNL, since it provides similar theoretical guarantees and has
almost the same global. FLECS uses sketching [38] and matrix compression for
efficient communication and allows for several Hessian learning techniques
(including truncated LSR1 update) and two iterate updates, namely Truncated
inverse Hessian approximation and FedSONIA, inspired by [30] and [22].

###  1.1 Contributions

We propose FLECS: a Federated LEarning Second-Order Framework via Compression
and Sketching with the following properties.

  * •

No additional storage on workers. In large-scale machine learning problems, it
is impossible to store the $d\times d$-dimensional Hessian approximations. For
a dataset with $d=50,000$ features, the local Hessian approximation requires
more than $10$ GB of memory. For FLECS, the Hessian approximations are stored
on the server, so memory requirements for the individual workers remain low.

  * •

Efficient communication via sketching and compression. To reduce the cost of
communication, we use a sketching technique. This allows the exchange of
$d\times m$ matrices instead of $d\times d$, where $m\ll d$ is the user
defined memory size, which could be set as $m=1$. Moreover, we exploit
compression to make the algorithm even more communication-efficient. This
strategy heavily improves the one in [36], since we decrease the dimension of
the transferred matrices $d/m$ times. Moreover, sketching reduces the number
of Hessian-vector products from $d$ to $m$.

  * •

Heterogeneous data setting. Our framework is flexible and does not assume that
data on the workers is independent and identically distributed.

  * •

Hessian Learning via Quasi-Newton updates. Motivated by quasi-Newton methods,
sketching techniques and inverse Hessian truncation, we propose the Truncated
L-SR1 update and Direct update for our Hessian approximations.

  * •

Theoretical Analysis. Convergence guarantees in both the strongly convex, and
nonconvex cases, are provided, as well as a local linear rate of convergence
under strong convexity independent of problem conditioning.

  * •

Competitive Numerical Experiments. Our numerical experiments highlight the
practical benefits of our proposed framework.

Table 1: Comparison between FedNL[36] and FLECS. | FedNL | FedNL-LS | FLECS
---|---|---|---
| [36] | [36] | [this paper]
heterogeneous data setting | ✓ | ✓ | ✓
additional storage on workers | $O(d^{2})$ | $O(d^{2})$ | $0$
supports biased and unbiased compression | ✓ | ✓ | ✓
communication cost (omitting compression) | $O(d^{2})$ | $O(d^{2})$ | $O(md)$
Hessian-vector products per iteration on worker | $d$ | $d$ | $m$
complexity per iteration on the server | $O(d^{3})$ | $O(d^{3})$ |  |  $O(d^{3})$ (Alg. 3)
---
$O(nmd^{2})$ (Alg. 4)
convergence in non convex case | ✗ | ✗ | ✓
strongly convex case: global convergence | ✗ | linear | linear
strongly convex case: local convergence | superlinear | ✗ | linear (Alg. 3)

#### Organization.

The remainder of the paper is organized as follows. Related works are
described in Section 2. Section 3 presents the proposed framework, including
all algorithmic details. In Section 4, we establish theoretical convergence
guarantees (all proofs can be found in the appendix). Numerical experiments
are provided in Section 5 (additional experiments are provided in the
appendix), and concluding remarks and potential future research directions are
given in Section 6.

##  2 Related Work

Second-order methods are well established in the optimization literature. Such
methods employ the Hessian (or an approximation to it) and thereby important
(possibly approximate and/or partial) curvature information is available for
utilization. The benefits include ‘better’ search directions, and faster rates
of convergence that may be independent of the condition number, although the
algorithms often have higher memory and computational requirements than first
order methods. To balance the trade off between improved convergence rates
versus increased costs, approaches that involve modifying or approximating the
Hessian have been proposed, and several of these are discussed now.
Arguably, some of the most popular algorithms that use approximations to the
Hessian are quasi-Newton methods. They employ information from the iterates to
build an approximation to the Hessian that improves in accuracy as the
algorithm progresses. Among these are algorithms that use BFGS, SR1, and DFP
updates [17, 29, 13], which are now considered to be classics in optimization
due to their effectiveness and practicality, and also because, under certain
conditions, theoretical results related to the accuracy of the Hessian
approximations are available. This motivates the current work, which
incorporates quasi-Newton type updates to learn a good approximation to the
Hessian.
Of particular interest is the SR1 update, which does not guarantee positive
definiteness of the Hessian approximation [13]. Previously, this was
considered to be a disadvantage, but in the case of non-convex optimization,
it would appear to be advantageous. L-SR1, a limited memory variant of SR1,
belongs to the class of limited memory quasi-Newton algorithms [25, 4]. These
methods store only the last $m$ gradient differences
$Y_{k}=[y_{k-m+1},\ldots,y_{k}]$ and iterate differences
$S_{k}=[s_{k-m+1},\ldots,s_{k}]$, where $m$ is the memory size, $y_{k}=\nabla
F(w_{k})-\nabla F(w_{k-1})$ and $s_{k}=w_{k}-w_{k-1}$. The compact form of the
L-SR1 update is [6]

$B_{k}=B_{0}+(Y_{k}-B_{0}S_{k})M_{k}^{-1}(Y_{k}-B_{0}S_{k})^{T},$ (2)

where $M_{k}=L_{k}+D_{k}+L_{k}^{T}-S_{k}^{T}B_{0}S_{k}$,
$~{}S_{k}^{T}Y_{k}=L_{k}+D_{k}+U_{k}$, $D_{k}$ is diagonal, $L_{k}$, $U_{k}$
are strictly lower and upper triangular, and $B_{0}$ is the initial
approximation.
The work [3] proposed the use of a sampled $S_{k}$, rather than one based upon
the previous $m$ iterate differences. The idea is, given the current iterate
$w_{k}$, simply sample $m$ random directions, and use the differences with
$w_{k}$ as the columns of $s_{k}$ (forming $Y_{k}$ appropriately). An $S_{k}$
with this form can be viewed as sketching matrix [38], and the Newton Sketch
algorithm in [31], shows that even randomly projected or a sub-sampled Hessian
provides enough curvature information for fast convergence.
The Nonconvex Newton method of [30] is a different approach, which is an
adaptation of Newton’s method for nonconvex problems. At each iteration, the
singular value decomposition of the Hessian is computed, small eigenvalues (in
absolute value) are truncated, and then a positive definite modification of
the Hessian is constructed. This is then used to create Newton-type direction.
SONIA [22] bridges the gap between first-order and second-order methods,
combining the ideas of quasi-Newton algorithms, sketching [31] and truncation
[30]. It constructs the search directions using curvature information in one
subspace and gradient information in the orthogonal complement.
Second-order methods for Federated Learning, can generally be divided in two
groups, based upon whether a homogeneous or heterogeneous data setting is
assumed. Algorithms in the first group rely on the concept of statistical
similarity, which, in practice, means that each local function $f_{i}$ is a
‘good’ approximation of the global objective function $F$. (More formally,
$\|\nabla^{2}F-\nabla^{2}f_{i}\|\leq\beta$ for some $\beta>0$.) Usually, in
this case, methods utilize the part of the data stored on the server and use
it to approximate the full Hessian [42, 15, 12, 5, 1].
No such statistical similarity (i.i.d.) assumption is made in the
heterogeneous setting, which is perhaps more realistic (e.g., it may be
impossible/impractical to have local data (on a mobile phone say), that
approximates the data on the server well). The main work in this setting is
FedNL [36]. The authors proposed a way to learn Hessian approximations and
then take an approximate Newton-type step. However, the main drawback is that
the Hessian approximations are stored locally and communications performed via
compressed $d\times d$ matrices. For large-scale optimization problems, it is
impossible to store the $d\times d$ matrix locally, because of memory
limitations. Extensions of FedNL can be found in [21], [32].
The approach proposed here combines ideas from quasi-Newton methods,
sketching, and Hessian learning, and develops a second-order Federated
Learning framework with efficient communications and without an additional
memory burden on workers. All Hessian approximations are stored on the server,
not the workers. Communications are based on a user-defined memory
hyperparameter $m$ and utilize compression to reduce communication. Sketching
allows to transfer $d\times m$ matrices and reduce the number of Hessian-
vector products on the worker. The proposed framework allows for two Hessian
learning techniques and two iterate update rules. The Truncated Inverse
Hessian approximation is based on the idea of truncating eigenvalues (as in
[30]) of Hessian approximate, while FedSONIA is based upon SONIA [22] and
reduces the iteration complexity of the step from $O(d^{3})$ to $O(nmd^{2})$.

##  3 FLECS

Inspired by the works [22, 36], we propose a general framework for second-
order methods for federated learning problems (1) that allows for different
forms for the iterates and different Hessian approximations/updates. One of
main goals of this framework is to create fast communication-efficient second-
order methods that do not place additional memory requirements on the workers.
To make broadcasting between the server and nodes cheap we employ sketching
and compression. Throughout the paper, $S_{k}\in\mathbb{R}^{d\times m}$
denotes a sketching matrix, where $m$ is the memory size ($m\ll d$). Any
operator that compresses a matrix is denoted by
$\mathcal{C}:\mathbb{R}^{d\times m}\to\mathbb{R}^{d\times m}$ , e.g. random
dithering, Rank-$R$, Top-$K$, Rand-$K$ (description of compressors can be
found in Appendix C). The Hessian approximation of the $i$-th worker at
iteration $k$ is denoted by $B_{k}^{i}$, and is stored on the central node.
Our framework can be described in three steps (further details follow in later
subsections):

  1. 1.

Worker step. Each worker receives the current iterate $w_{k}$, and their local
Hessian approximation multiplied by a sketching matrix, i.e., the product
$B_{k}^{i}S_{k}\in\mathbb{R}^{d\times n}$. The matrix is sent directly (i.e.,
no compression), since, typically, downloading is faster than uploading. Each
worker then computes their local Hessian-sketch product
$Y_{k}^{i}:=\nabla^{2}f_{i}(w_{k})S_{k}$ and sends back the compressed
difference
$C_{k}^{i}:=\mathcal{C}(Y_{k}^{i}-B_{k}^{i}S_{k})\in\mathbb{R}^{d\times m}$.

  2. 2.

Hessian approximation update. The server receives $C_{k}^{i}$,
$M_{k}^{i}:=S_{k}^{T}(\nabla^{2}f_{i}(w_{k}))S_{k}$, and $\nabla f_{i}(w_{k})$
from all $i=1,\ldots,n$ workers. Using this information, as well as the
current (local) Hessian approximation $B_{k}^{i}$, the new approximation
$B_{k+1}^{i}$ is formed. Our framework allows $B_{k+1}^{i}$ to be computed by
a (1) Truncated L-SR1 update, or (2) Direct update.

  3. 3.

Iterate update. The server forms the matrices

$B_{k}:=\tfrac{1}{n}\textstyle{\sum}_{i=1}^{n}B_{k}^{i},~{}~{}\nabla
F(w_{k}),~{}~{}M_{k}=\frac{1}{n}\textstyle{\sum}_{i=1}^{n}M_{k}^{i},$

and calculates the new iterate $w_{k+1}$ via update rule:

$w_{k+1}=w_{k}+\alpha_{k}p_{k}$ (3)

where $\alpha_{k}>0$ is the step-size. There are several options for computing
the search direction $p_{k}$, including (1) a truncated inverse Hessian
approximation step, or (2) via a SONIA step for Federated Learning (FedSONIA).
The process then repeats from Step 1.

###  3.1 Initialization

The algorithm is initialized with a user defined memory size $m\ll d$, and $n$
matrices $B_{0}^{i}\in\mathbb{R}^{d\times d}$ are stored on the central node,
each representing an approximation to the local Hessian for the $i$th worker.
Various additional parameters are initialized depending on the selected update
rule (to be described later). At iteration $k>0$ of FLECS, $m$ directions
$\left\\{s_{1},\ldots,s_{m}\right\\}$ are randomly sampled at all nodes
(including the server) and combined together, forming the columns of a
sketching matrix $S_{k}\in\mathbb{R}^{d\times m}$. Note that $S_{k}$ is the
same for all machines and the main node; this is _guaranteed_ by setting the
random seed to be equal to the iteration number $k$ during sampling.

###  3.2 Worker Step

Each of the workers receives $B_{k}^{i}S_{k}$, and the iterate $w_{k}$, from
the server. On each node we compute

$C_{k}^{i}:=\mathcal{C}(Y_{k}^{i}-B_{k}^{i}S_{k}),~{}\nabla
f_{i}(w_{k}),~{}~{}M_{k}^{i}:=S_{k}^{T}Y_{k}^{i}=S_{k}^{T}(\nabla^{2}f_{i}(w_{k}))S_{k},~{}~{}i=1,\ldots,n.$
(4)

The quantities in (4) are used to construct the Hessian approximation on the
main node. Since $M_{k}^{i}$ is an $m\times m$, it can be sent directly
without compression. So too can the vector $\nabla f_{i}(w_{k})$. For the
matrix $Y_{k}^{i}$ we use the compression operator $\mathcal{C}$. In order to
make the variance of this compression low, we use an error-feedback technique
and send $C_{k}^{i}=\mathcal{C}(Y_{k}^{i}-B_{k}^{i}S_{k}).$

###  3.3 Hessian Approximation update

The server receives (4) from each worker $i=1,\ldots,n$. Firstly, we restore
an approximation to $Y_{k}^{i}$ as

$\widetilde{Y}_{k}^{i}=C_{k}^{i}+B_{k}^{i}S_{k}.$ (5)

Then, we apply one of the following options to update the Hessian
approximation.

Hessian update option 1. Truncated L-SR1 update (Algorithm 1). In this case we
compute the approximate Hessian of the $i$-th worker using the L-SR1 update:

$B_{k+1}^{i}=B_{k}^{i}+(\widetilde{Y}_{k}-B_{k}^{i}S_{k})(M_{k}^{i}-S_{k}^{T}B_{k}^{i}S_{k})^{\dagger}(\widetilde{Y}_{k}-B_{k}^{i}S_{k})^{T}.$
(6)

To ensure numerically stability of the update (6), truncation is employed. At
first, we take the spectral decomposition of the small $(m\times m)$ matrix

$M_{k}^{i}-S_{k}^{T}B_{k}^{i}S_{k}=U_{k}^{i}L_{k}^{i}(U_{k}^{i})^{T},$

where $L_{k}^{i}\in\mathbb{R}^{d\times d}$ is a diagonal matrix containing the
eigenvalues, and $U_{k}^{i}$ is an orthogonal matrix containing the
corresponding eigenvectors. Next, we compute the inverse $(L_{k}^{i})^{-1}$,
and truncate small values. In particular, define
$(l_{k}^{i})_{j}^{-1}:=(L_{k}^{i})^{-1}_{jj}$, and set all elements such that
$|(l_{k}^{i})_{j}^{-1}|\leq\omega$ to be $0$, where $\omega>0$ is a user
defined parameter that must be initialized. The truncated matrix is denoted by
$[(L_{k}^{i})^{-1}]_{\omega}$. Finally, the update (6) becomes

$B_{k+1}=B_{k}+(\widetilde{Y}_{k}-B_{k}^{i}S_{k})U_{k}^{i}[(L_{k}^{i})^{-1}]_{\omega}(U_{k}^{i})^{T}(\widetilde{Y}_{k}-B_{k}^{i}S_{k})^{T}.$
(7)

Algorithm 1 Truncated L-SR1 update

0: $\widetilde{Y}_{k}^{i}\in\mathbb{R}^{d\times
m},~{}M_{k}^{i}\in\mathbb{R}^{m\times m},~{}B_{k}^{i}\in\mathbb{R}^{d\times
d},~{}S_{k}\in\mathbb{R}^{d\times m}$ for $i=1,\ldots,n$, $~{}\omega>0$ –
truncation constant.

1: On the server:

2: for $i=1,\ldots,n$ do

3: compute
$(M_{k}^{i}-(S_{k}^{i})^{T}\widetilde{Y}_{k}^{i})=U_{k}^{i}L_{k}^{i}(U_{k}^{i})^{T}$;

4: truncate $(L_{k}^{i})^{-1}$ to form $[(L_{k}^{i})^{-1}]_{\omega}$;

5: compute $B_{k+1}^{i}$ via (7).

6: end for

Hessian update option 2. Direct update (Algorithm 2). Setting $B_{0}=0$ in (6)
gives

$\widetilde{B}^{i}_{k}=\widetilde{Y}^{i}_{k}(M^{i}_{k})^{\dagger}(\widetilde{Y}^{i}_{k})^{T},$
(8)

which we refer to as the Direct update. Then, we utilize the following
learning technique

$B_{k+1}^{i}=(1-\beta_{k})B_{k}^{i}+\beta_{k}\widetilde{B}^{i}_{k},$ (9)

where $0<\beta_{k}\leq 1$ is the learning rate. Thus, for the Direct update,
the new Hessian approximation is a convex combination of the previous update
and the matrix in (8).
Note, that for both Option 1 and Option 2, the Hessian approximation
$B_{k+1}^{i}$ is symmetric, since $M_{k}^{i}=S_{k}^{T}\nabla
f_{i}(w_{k})S_{k}$, $B_{k}^{i}$, and $M_{k}^{i}-S_{k}^{T}B_{k}^{i}S_{k}$ are
symmetric for $i=1,\ldots,n$.

Algorithm 2 Direct update.

0: $\widetilde{Y}_{k}^{i}\in\mathbb{R}^{d\times
m},~{}M_{k}^{i}\in\mathbb{R}^{m\times m},~{}B_{k}^{i}\in\mathbb{R}^{d\times
d}$ $\forall i$

1: On the server: $0<\beta_{k}\leq 1$ – learning rate.

2: for $i=1,\ldots,n$ do

3: compute
$\widetilde{B}^{i}_{k}=\widetilde{Y}^{i}_{k}(M^{i}_{k})^{\dagger}(\widetilde{Y}^{i}_{k})^{T}$;

4: select learning rate $\beta_{k}$

5: compute
$B^{i}_{k+1}=(1-\beta_{k})B^{i}_{k}+\beta_{k}\widetilde{B}^{i}_{k}$.

6: end for

###  3.4 Iterate update

Note that, at the start of the iterate update step, all the Hessian
approximations $B_{k+1}^{i}$, and all variables from (4), are available on the
main node. Then, the server forms

$\displaystyle\nabla F(w_{k})$
$\displaystyle=\tfrac{1}{n}\textstyle{\sum}_{i=1}^{n}f_{i}(w_{k}),$
$\displaystyle\widetilde{Y}_{k}$
$\displaystyle:=\tfrac{1}{n}\textstyle{\sum}\widetilde{Y}_{k}^{i}=\tfrac{1}{n}\textstyle{\sum}_{i=1}^{n}(C_{k}^{i}+H_{k}^{i}S_{k}),$
$\displaystyle M_{k}$
$\displaystyle:=\tfrac{1}{n}\textstyle{\sum}_{i=1}^{n}M_{k}^{i}=S_{k}\nabla^{2}F(w_{k})S_{k}^{T},$
$\displaystyle B_{k+1}$
$\displaystyle:=\tfrac{1}{n}\textstyle{\sum}_{i=1}^{n}B_{k+1}^{i}.$

Now the goal is to compute a search direction $p_{k}$, of the general form:

$p_{k}=-A_{k}\nabla F(w_{k}),$ (10)

for some matrix $A_{k}$. Two specific options for $A_{k}$ are described now.

Search direction option 1. Truncated inverse Hessian approximation (Algorithm
3). Here, $A_{k}$ is taken to be the truncated inverse Hessian approximation
$\left(|B_{k+1}|_{\omega}^{\Omega}\right)^{-1}$. Consider the following.

######  Definition 1.

Let $B_{k},V_{k},\Lambda_{k}$ be matrices such that
$B_{k}=V_{k}\Lambda_{k}V_{k}^{T}$, and let $0<\omega\leq\Omega$. The truncated
inverse Hessian approximation of $B_{k}$ is
$\left(|B_{k}|_{\omega}^{\Omega}\right)^{-1}:=V_{k}(|\Lambda_{k}|_{\omega}^{\Omega})^{-1}V_{k}^{T}$,
where
$(|\Lambda_{k}|_{\omega}^{\Omega})_{ii}=\min\left\\{\max\left\\{|\Lambda|_{ii},\omega\right\\},\Omega\right\\}.$

Definition 1 was proposed in [30] and was used to provide a convergence
guarantee for their Nonconvex Newton method (to a local minimum). Firstly, an
eigen-decomposition of $B_{k}$ is computed, but with every eigenvalue replaced
by its absolute value. Secondly, a thresholding step is applied, so that any
eigenvalue (in absolute value) that is smaller (resp. greater) than a user
defined threshold $\omega$ (resp. $\Omega$) is replaced by $\omega$ (resp.
$\Omega$). This truncation is crucial for several reasons. Firstly, it ensures
$B_{k}$ is well-defined, which leads to stable updates. Secondly, it ensures
that after truncation, the resulting matrix will be full-rank and positive
definite. Therefore, it guarantees that the following search direction is a
_decent_ direction

$p_{k}=\left(|B_{k+1}|_{\omega}^{\Omega}\right)^{-1}\nabla F(w_{k}).$ (11)

Algorithm 3 Truncated inverse Hessian approximation

0: $\nabla
F(w_{k})\in\mathbb{R}^{d},~{}\widetilde{Y}_{k}\in\mathbb{R}^{d\times
m},~{}M_{k}\in\mathbb{R}^{m\times m},~{}B_{k+1}\in\mathbb{R}^{d\times d}$,
$\Omega>\omega>0$ – truncation constants.

1: On the server:

2: compute spectral decomposition $B_{k+1}=V_{k}\Lambda_{k}V_{k}^{T}$;

3: truncate $\Lambda_{k}$ to form $|\Lambda_{k}|_{\omega}^{\Omega}$ via
Definition 1;

4: compute search direction $p_{k}$ via (11);

5: return $p_{k}$.

Search direction option 2. FedSONIA (Algorithm 4). FedSONIA is a modification
of SONIA [22] that is adapted to the Federated Learning setting. The key idea
of SONIA is to utilize curvature information in one subspace and perform a
gradient step in orthogonal complement. The intuition is that, by
incorporating some curvature information the resulting search direction should
be ‘better’ than a pure gradient direction, but the user has control over the
dimension of each subspace, and therefore they have control over the cost of
obtaining the curvature information. Truncation (Definition 1) and modified
SR1 update (eq. (6)) ensure that the Hessian-approximation is full-rank,
positive and well defined, which is crucial to avoid checking condition on the
curvature pairs $Y_{k}^{T}S_{k}$ ($\widetilde{Y}_{k}^{T}S_{k}$ for FedSONIA).
To construct the approximation of the Hessian we use (8). Note, that SONIA
does not need the Hessian approximation $B_{k}$ itself to commit a step.
Instead, it utilizes only Hessian-matrix product $Y_{k}$ and sketching matrix
$S_{k}$.
Define
$B_{k+1}^{\text{SONIA}}=\widetilde{Y}_{k}(M_{k})^{\dagger}\widetilde{Y}_{k}^{T}$.
Following the strategy in [22], apply an economy $QR$ factorization of
$\widetilde{Y}_{k}=Q_{k}R_{k}$, so that $Q_{k}\in\mathbb{R}^{d\times m}$ has
orthonormal columns and $R_{k}\in\mathbb{R}^{m\times m}$ is upper triangular.
Thus, $B_{k+1}^{\text{SONIA}}=Q_{k}R_{k}(M_{k})^{\dagger}R_{k}^{T}Q_{k}^{T}$.
Applying the spectral decomposition, we obtain
$R_{k}(M_{k})^{\dagger}R_{k}^{T}=V_{k}\Lambda_{k}V_{k}^{T}$, where the columns
of $V_{k}\in\mathbb{R}^{m\times m}$ form an orthonormal basis (of
eigenvectors), and $\Lambda_{k}\in\mathbb{R}^{m\times m}$ is a diagonal matrix
containing the corresponding eigenvalues. Therefore,

$B_{k+1}^{\text{SONIA}}=Q_{k}V_{k}\Lambda_{k}V_{k}^{T}Q_{k}^{T}\stackrel{{\scriptstyle\widetilde{V}_{k}:=Q_{k}V_{k}}}{{=}}\widetilde{V}_{k}\Lambda_{k}\widetilde{V}_{k}^{T},$

where $\widetilde{V}_{k}\in\mathbb{R}^{d\times m}$ has orthonormal columns
since $Q_{k}$ has orthonormal columns and $V_{k}$ is an orthogonal matrix.
Then, we apply Definition 1 to form the truncated approximation

$\left(|B_{k+1}^{\text{SONIA}}|_{\omega}^{\Omega}\right)^{-1}=\widetilde{V}_{k}(|\Lambda_{k}|_{\omega}^{\Omega})^{-1}\widetilde{V}_{k}^{T}.$
(12)

Next, we perform subspace decomposition on the current gradient direction:

$\nabla F(w_{k})=g_{k}+g_{k}^{\perp},$ (13)

where $g_{k}=\widetilde{V}_{k}\widetilde{V}_{k}^{T}\nabla F(w_{k})$ and
$g_{k}^{\perp}=(I-\widetilde{V}_{k}\widetilde{V}_{k}^{T})\nabla F(w_{k})$.
Putting this all together, the FedSONIA search direction (similar to [22]) is:

$p_{k}:=-\left(|B_{k+1}^{\text{SONIA}}|_{\omega}^{\Omega}\right)^{-1}g_{k}-\rho_{k}g^{\perp}_{k}\stackrel{{\scriptstyle\eqref{eq:A_sonia},~{}\widetilde{V}_{k}^{T}V_{k}=1}}{{=}}-\left(\left(|B_{k+1}^{\text{SONIA}}|_{\omega}^{\Omega}\right)^{-1}+\rho_{k}(I-\widetilde{V}_{k}\widetilde{V}_{k}^{T})\right)\nabla
F(w_{k}),$ (14)

where
$\rho_{k}\in[1/\Omega,1/(\max\limits_{i}[(|\Lambda_{k}|_{\omega}^{\Omega})^{-1}]_{ii})]$.
One has control over the cost of the curvature information in (14) via the
memory size $m$. Note that this will also affect the communication complexity.

Algorithm 4 FedSONIA

0: $\nabla
F(w_{k})\in\mathbb{R}^{d},~{}\widetilde{Y}_{k}\in\mathbb{R}^{d\times
m},~{}M_{k}\in\mathbb{R}^{m\times m}$, $\Omega>\omega>0$ – truncation
constants.

1: On the server:

2: compute
$B_{k}^{\text{SONIA}}:=\widetilde{Y}_{k}(M_{k})^{\dagger}\widetilde{Y}_{k}^{T}$;

3: compute $QR$ factorization of $\widetilde{Y}_{k}(=Q_{k}R_{k})$;

4: compute spectral decomposition of
$R_{k}(M_{k})^{\dagger}R_{k}^{T}(=V_{k}\Lambda_{k}V_{k}^{T})$;

5: construct $\widetilde{V}_{k}:=Q_{k}V_{k}$;

6: truncate $\Lambda_{k}$ to form $|\Lambda_{k}|_{\omega}^{\Omega}$ via
Definition 1;

7: Set $\rho_{k}$ and decompose gradient via (13);

8: Compute search direction $p_{k}$ via (14);

9: return $p_{k}$.

Finally, a new search direction $p_{k}$ is computed, which depends on the
chosen update strategy (Algorithm 3 or 4) and to compute the new iterate
$w_{k+1}$ via (3). Recall that for both these strategies the resulting search
direction $p_{k}$ has the form (10), where
$A_{k}=(|B_{k+1}|_{\omega}^{\Omega})^{-1}$ for the Truncated inverse Hessian
approximation (11), and
$A_{k}=\left(|B_{k+1}^{\text{SONIA}}|_{\omega}^{\Omega}\right)^{-1}+\rho_{k}(I-\widetilde{V}_{k}\widetilde{V}_{k}^{T})$
for FedSONIA (14). The process then restarts (from Section 3.2). Our framework
is summarized in Algorithm 5.
Complexity. The per iteration complexity of FLECS depends on the chosen
options for both the Hessian approximation and the iterate updates. The
complexity of the Truncated L-SR1 update for one Hessian approximation
consists of matrix products $(O(md^{2}))$ and a singular value decomposition
of $m\times m$ matrix $(O(m^{3}))$. For the Direct update, the matrix products
cost $O(md^{2})$ and forming the pseudo-inverse of the $m\times m$ matrix
costs $O(m^{3})$. Therefore, the total complexity of either Hessian
approximation update is $O(nmd^{2})$. The complexity of computing search
direction via Truncated inverse Hessian approximation is $O(d^{3})$, since we
utilize spectral decomposition and products of $d\times d$ matrices. FedSONIA
update comprises 1) QR factorization $(O(dm^{2}))$; 2) pseudoinverse of
$m\times m$ matrix $(O(m^{3}))$; 3) spectral decomposition of $m\times m$
matrix $(O(m^{2}))$; which gives total complexity $O(dm^{2})$. Finally, the
complexity of FLECS on the server only depends on options for the search
direction and equals $O(d^{3}+nmd^{2})$ for Truncated Inverse Hessian
approximation and $O(nmd^{2})$ for FedSONIA. The the worker step’s complexity
consists of $m$ Hessian-vector products and matrix multiplication
$(O(md^{2}))$.

Algorithm 5 FLECS

0: $w_{0}$ – starting point, $m$ – memory size, $B_{0}^{i}$ – initial Hessian
approximations for each worker $i=1\ldots n$ on the server, $0<\omega<\Omega$
– truncation constants.

1: for $k=0,1,\ldots$ do

2: On $i$-th machine:

3: sample $S_{k}\in\mathbb{R}^{d\times m}$;

4: collect $H_{k}^{i}S_{k},~{}w_{k}$ from the server;

5: compute
$Y_{k}^{i}:=\nabla^{2}f_{i}(w_{k})S_{k},~{}U_{k}^{i}:=S_{k}^{T}Y_{k}^{i}$;

6: send $\nabla
f_{i}(w_{k}),~{}M_{k}^{i},~{}C_{k}^{i}=\mathcal{C}(Y_{k}^{i}-H_{k}^{i}S_{k})$
to the server.

7: On the server:

8: sample $S_{k}$;

9: collect $C_{k}^{i},~{}M_{k}^{i},~{}i=1\ldots n$ from workers;

10: compute $\widetilde{Y}_{k}^{i}=C_{k}^{i}+H_{k}^{i}S_{k}$;

11: compute $B_{k+1}^{i}$ via Algorithm 1 or select learning rate $\beta_{k}$
and compute $B_{k+1}^{i}$ via Algorithm 2;

12: form $\nabla F(w_{k}),\widetilde{Y}_{k},M_{k},B_{k+1}$ as average over
workers of $\nabla
f_{i}(w_{k}),\widetilde{Y}_{k}^{i},M^{i}_{k},B^{i}_{k+1},~{}i=1,\ldots,n$;

13: compute search direction $p_{k}$ via Algorithm 3 or 4;

14: select stepsize $\alpha_{k}$ and set $w_{k+1}=w_{k}+\alpha_{k}p_{k}$;

15: sample $S_{k+1}\in\mathbb{R}^{d\times m}$;

16: send $w_{k},B_{k+1}S_{k+1}$ to all workers.

17: end for

##  4 FLECS: Convergence Theory

In this section, convergence theory for our FLECS framework is presented. We
provide global convergence results in both the strongly convex, and nonconvex
cases. A local convergence result is also presented, which shows local linear
convergence of the iterates under strong convexity. We define $w_{0}$ to be
the initial point, $w^{\star}=\arg\min\limits_{w\in\mathbb{R}^{d}}F(w)$ and
$F^{\star}=F(w^{\star})$.

######  Assumption 1.

The function $F(w)$ is twice continuously differentiable for all
$w\in\mathbb{R}^{d}$.

We begin by listing the preliminary lemma.

######  Lemma 1.

Let Assumption 1 hold. Let $A_{k}$ be defined as in (10) depending on the
choice of iterates update. Then, for any $k\geq 0$ there exist constants
$0<\mu_{1}\leq\mu_{2}$ such that $\mu_{1}I\preceq A_{k}\preceq\mu_{2}I$.

We introduce assumptions and theoretical results for the strongly convex case.

######  Assumption 2.

There exist positive constants $\mu$ and $L$ such that $\mu
I\preceq\nabla^{2}F(w)\preceq LI,$ $\forall w\in\mathbb{R}^{d}$.

The following theorem establishes global linear convergence of FLECS under
strong convexity.

######  Theorem 1.

Let Assumptions 1 and 2 hold. Then for iterates $\\{w_{k}\\}$ generated by
Algorithm 5 with $0<\alpha_{k}=\alpha\leq\frac{\mu_{1}}{\mu_{2}^{2}L}$ we have
$F(w_{k})-F^{\star}\leq(1-\alpha\mu\mu_{1})^{k}[F(w_{0})-F^{\star}]~{}$ for
all $k\geq 0$.

We are now ready to present a convergence guarantee which shows local linear
convergence for FLECS under strong convexity. Note that, while Theorem 1 shows
that FLECS converges linearly under strong convexity, the rate depends on
parameters related to the Hessian approximations. The following local rate is
_independent_ of problem conditioning.

######  Theorem 2.

Let Assumptions 1 and 2 hold, and $\|w^{0}-w^{*}\|\leq\frac{\mu^{2}}{2}$,
$\frac{1}{n}\sum\limits_{i=1}^{n}\|B_{k+1}^{i}-\nabla^{2}f_{i}(w^{*})\|_{\text{F}}\leq\frac{2\mu^{2}}{L^{2}}$.
Then for iterates $\\{w_{k}\\}$ generated by Algorithm 5 with Truncated
inverse Hessian approximation step (Algorithm 3) with parameters
$\alpha_{k}=\alpha=1,~{}0\leq\omega\leq\mu,~{}\Omega\geq L$ we have

$\|w_{k}-w^{*}\|^{2}\leq\tfrac{1}{2^{k}}\|w_{0}-w^{*}\|^{2}$

for all $k\geq 0$.

We introduce the assumptions for nonconvex case.

######  Assumption 3.

The function $F$ is bounded below by a scalar $\hat{F}$.

######  Assumption 4.

The gradients of $F$ are $L$-Lipschitz continuous for all
$w\in\mathbb{R}^{n}$.

The following result shows that FLECS is convergent when problem (1) is
nonconvex.

######  Theorem 3.

Let Assumptions1, 3, 4 hold and $w_{0}$ be the starting point. Then after
$T>0$ iterations of Algorithm 5 with
$0<\alpha_{k}=\alpha\leq\frac{\mu_{1}}{\mu_{2}^{2}L}$, we have

$\tfrac{1}{T}\textstyle{\sum}_{k=0}^{T-1}\|\nabla
F(w_{k})\|^{2}\leq\tfrac{2[F(w_{0})-\hat{F}]}{\alpha\mu_{1}T}\stackrel{{\scriptstyle
T\to\infty}}{{\longrightarrow}}0.$

##  5 Numerical Experiments

Figure 1: Comparison of objective function $F(w_{k})$ and the squared norm of
gradient $\|\nabla F(w_{k})\|^{2}$.

We present numerical experiments for FLECS on a regularized logistic
regression problem:

$\min\left\\{\tfrac{1}{n}\textstyle{\sum}_{i=1}^{n}\tfrac{1}{r}\textstyle{\sum}_{j=1}^{r}\log(1+\exp(-b_{ij}a_{ij}^{T}w))+\mu\|w\|^{2}\right\\},$

where $\left\\{a_{ij},b_{ij}\right\\}_{j\in[m]}$ are data points on $i$-th
device. We use three datasets from the LIBSVM library [7]: a9a ($123$
features), gisette-scale ($5000$ features) and real-sim ($20958$ features). We
use Top-$K$ with $k=4d$ and random dithering [2] applied to each matrix column
with $s=128$ levels and $\infty$-norm as a matrix compressors. The performance
of FLECS is compared with FedNL.
In our experiments we set regularization parameter $\mu=10^{-5}$, the memory
size to $m=16$ and use $B_{k}^{i}=0$ as the initial Hessian approximation for
FLECS and FedNL for random dithering and $B_{k}^{i}=\nabla f_{i}(w_{0})$ for
Top-$K$. Truncation parameters were set to $\omega=10^{-3},~{}\Omega=10^{8}$,
and for the FedSONIA update we set $\rho=\frac{1}{\Omega}$ and $\alpha=1$. For
FLECS we use a Truncated L-SR1 Update (Algorithm 1) and a Direct update
(Algorithm 2) with learning rate $\beta=1$ as the Hessian learning techniques.
For FedNL we use the same $\beta$.
The following three rules were used to compute the search direction $p_{k}$:
(1) a Truncated Inverse Hessian Approximation (Algorithm 3), (2) FedSONIA
(Algorithm 4) with $\rho_{k}=\frac{1}{\Omega}$, and (3) a regularized update
where $B_{k+1}=B_{k}+\frac{1}{n}\sum_{i=1}^{n}\|C_{k}^{i}\|I$, where
$C_{k}^{i}$ are defined in (4). The last update is motivated by [36]. Note,
that in the case of FLECS this update does not guarantee positive definiteness
of Hessian approximation $B_{k}$ if $m<d$, but it can still perform well if
$m$ is large enough.
For FedNL after random dithering the matrix $B_{k+1}$ might not be symmetric,
so we use $\frac{1}{2}(B_{k+1}+B_{k+1}^{T})$ as Hessian approximation. In this
case, we compare algorithms on the number of bits sent to the server per
iteration. For Top-K compression, the communication complexity is roughly the
same, so we estimate the performance based on the number of Hessian-vector
products.
For the iterates update (3) we set the step-size $\alpha=1$ for all algorithms
and datasets except FedSONIA on gisette and real-sim, where we set
$\alpha=0.1$, as it gives the best performance. Further details specifying
algorithm parameters and implementation, as well as additional numerical
experiments, can be found in Appendix B.

##  6 Conclusion

This work presents a new second order framework for Federated Learning
problems. The framework employs sketching and compression to ensure
communication costs are low, and Hessian approximations are stored on the
central server so that memory requirements for the workers is low. Due to the
inclusion of partial approximate second order information, the algorithm
enjoys favourable convergence guarantees, and numerical experiments show the
practical benefits of our approach.

## References

  * [1] A. Agafonov, P. Dvurechensky, G. Scutari, A. Gasnikov, D. Kamzolov, A. Lukashevich, and A. Daneshmand.  An accelerated second-order method for distributed stochastic optimization.  In 2021 60th IEEE Conference on Decision and Control (CDC), pages 2407–2413, 2021.
  * [2] D. Alistarh, D. Grubic, J. Li, R. Tomioka, and M. Vojnovic.  Qsgd: Communication-efficient sgd via gradient quantization and encoding.  Advances in Neural Information Processing Systems, 30, 2017.
  * [3] A. S. Berahas, M. Jahani, P. Richtárik, and M. Takáč.  Quasi-newton methods for machine learning: forget the past, just sample.  Optimization Methods and Software, pages 1–37, 2021.
  * [4] J. Brust, J. B. Erway, and R. F. Marcia.  On solving l-sr1 trust-region subproblems.  Computational Optimization and Applications, 66(2):245–266, 2017\.
  * [5] B. Bullins, K. Patel, O. Shamir, N. Srebro, and B. E. Woodworth.  A stochastic newton algorithm for distributed convex optimization.  Advances in Neural Information Processing Systems, 34, 2021.
  * [6] R. H. Byrd, J. Nocedal, and R. B. Schnabel.  Representations of quasi-newton matrices and their use in limited memory methods.  Mathematical Programming, 63(1):129–156, 1994.
  * [7] C.-C. Chang and C.-J. Lin.  Libsvm: a library for support vector machines.  ACM transactions on intelligent systems and technology (TIST), 2(3):1–27, 2011.
  * [8] D. Chen, C. S. Hong, Y. Zha, Y. Zhang, X. Liu, and Z. Han.  Fedsvrg based communication efficient scheme for federated learning in mec networks.  IEEE Transactions on Vehicular Technology, 70(7):7300–7304, 2021\.
  * [9] Y. Chen, R. S. Blum, M. Takac, and B. M. Sadler.  Distributed learning with sparsified gradient differences.  IEEE Journal of Selected Topics in Signal Processing, 2022.
  * [10] J. K. cný, H. B. McMahan, D. Ramage, and P. Richtarik.  Federated optimization: Distributed machine learning for on-device intelligence.  arXiv preprint arXiv:1610.02527, 2016.
  * [11] J. K. cný, H. B. McMahan, F. X. Yu, P. Richtarik, A. T. Suresh, and D. Bacon.  Federated learning: Strategies for improving communication efficiency.  In NIPS Workshop on Private Multi-Party Machine Learning, 2016.
  * [12] A. Daneshmand, G. Scutari, P. Dvurechensky, and A. Gasnikov.  Newton method over networks is fast up to the statistical precision, 2021\.
  * [13] J. E. Dennis, Jr and J. J. Moré.  Quasi-newton methods, motivation and theory.  SIAM review, 19(1):46–89, 1977.
  * [14] P. Dvurechenskii, D. Dvinskikh, A. Gasnikov, C. Uribe, and A. Nedich.  Decentralize and randomize: Faster algorithm for wasserstein barycenters.  In Advances in Neural Information Processing Systems, pages 10760–10770, 2018.
  * [15] P. Dvurechensky, D. Kamzolov, A. Lukashevich, S. Lee, E. Ordentlich, C. A. Uribe, and A. Gasnikov.  Hyperfast second-order local solvers for efficient statistically preconditioned distributed optimization.  arXiv preprint arXiv:2102.08246, 2021.
  * [16] P. K. et al.  Advances and open problems in federated learning, 2019.
  * [17] R. Fletcher.  Practical methods of optimization.  John Wiley & Sons, 2013.
  * [18] H. He, X. Yu, J. Zhang, S. Song, and K. B. Letaief.  Cell-free massive mimo for 6g wireless communication networks.  Journal of Communications and Information Networks, 6(4):321–335, 2021.
  * [19] S. Horváth, C.-Y. Ho, L. Horvath, A. N. Sahu, M. Canini, and P. Richtárik.  Natural compression for distributed deep learning.  arXiv preprint arXiv:1905.10988, 2019.
  * [20] S. Horváth, D. Kovalev, K. Mishchenko, S. Stich, and P. Richtárik.  Stochastic distributed learning with gradient quantization and variance reduction, 2019.
  * [21] R. Islamov, X. Qian, and P. Richtárik.  Distributed second order methods with fast rates and compressed communication.  In International Conference on Machine Learning, pages 4617–4628. PMLR, 2021.
  * [22] M. Jahani, M. Nazari, R. Tappenden, A. Berahas, and M. Takác.  Sonia: A symmetric blockwise truncated optimization algorithm.  In International Conference on Artificial Intelligence and Statistics, pages 487–495. PMLR, 2021.
  * [23] T. Kraska, A. Talwalkar, J. C. Duchi, R. Griffith, M. J. Franklin, and M. I. Jordan.  Mlbase: A distributed machine-learning system.  In CIDR, volume 1, pages 2–1, 2013.
  * [24] M. Li, D. G. Andersen, J. W. Park, A. J. Smola, A. Ahmed, V. Josifovski, J. Long, E. J. Shekita, and B.-Y. Su.  Scaling distributed machine learning with the parameter server.  In 11th USENIX Symposium on Operating Systems Design and Implementation (OSDI 14), pages 583–598, 2014.
  * [25] X. Lu.  A study of the limited memory SR1 method in practice.  University of Colorado at Boulder, 1996.
  * [26] H. B. McMahan, E. Moore, D. Ramage, S. Hampson, and B. A. y Arcas.  Communication-efficient learning of deep networks from decentralized data.  In In Proceedings of the 20th International Conference on Artificial Intelligence and Statistics (AISTATS), 2017.
  * [27] K. Mishchenko, E. Gorbunov, M. Takáč, , and P. Richtárik.  Distributed learning with compressed gradient differences, 2019.
  * [28] A. Nedić, A. Olshevsky, and C. A. Uribe.  Fast convergence rates for distributed non-bayesian learning.  IEEE Transactions on Automatic Control, 62(11):5538–5553, 2017\.
  * [29] J. Nocedal and S. J. Wright.  Numerical optimization.  Springer, 1999.
  * [30] S. Paternain, A. Mokhtari, and A. Ribeiro.  A newton-based method for nonconvex optimization with fast evasion of saddle points.  SIAM Journal on Optimization, 29(1):343–368, 2019.
  * [31] M. Pilanci and M. J. Wainwright.  Newton sketch: A near linear-time optimization algorithm with linear-quadratic convergence.  SIAM Journal on Optimization, 27(1):205–245, 2017.
  * [32] X. Qian, R. Islamov, M. Safaryan, and P. Richtárik.  Basis matters: Better communication-efficient second order methods for federated learning.  arXiv preprint arXiv:2111.01847, 2021.
  * [33] M. Rabbat and R. Nowak.  Decentralized source localization and tracking wireless sensor networks.  In Proceedings of the IEEE International Conference on Acoustics, Speech, and Signal Processing, volume 3, pages 921–924, 2004.
  * [34] V. Rajaravivarma, Y. Yang, and T. Yang.  An overview of wireless sensor network and applications.  In Proceedings of the 35th Southeastern Symposium on System Theory, 2003., pages 432–436. IEEE, 2003.
  * [35] S. S. Ram, V. V. Veeravalli, and A. Nedic.  Distributed non-autonomous power control through distributed convex optimization.  In IEEE INFOCOM 2009, pages 3001–3005. IEEE, 2009.
  * [36] M. Safaryan, R. Islamov, X. Qian, and P. Richtárik.  Fednl: Making newton-type methods applicable to federated learning.  arXiv preprint arXiv:2106.02969, 2021.
  * [37] C. A. Uribe, D. Dvinskikh, P. Dvurechensky, A. Gasnikov, and A. Nedić.  Distributed computation of Wasserstein barycenters over networks.  In 2018 IEEE 57th Annual Conference on Decision and Control (CDC), 2018.  Accepted, arXiv:1803.02933.
  * [38] D. P. Woodruff.  Sketching as a tool for numerical linear algebra.  arXiv preprint arXiv:1411.4357, 2014.
  * [39] L. Xiao and S. Boyd.  Optimal scaling of a gradient method for distributed resource allocation.  Journal of Optimization Theory and Applications, 129(3):469–488, 2006.
  * [40] C. Xie, O. Koyejo, I. Gupta, and H. Lin.  Local adaalter: Communicationefficient stochastic gradient descent with adaptive learning rates, 2019.
  * [41] J. Yick, B. Mukherjee, and D. Ghosal.  Wireless sensor network survey.  Computer networks, 52(12):2292–2330, 2008.
  * [42] Y. Zhang and L. Xiao.  Communication-Efficient Distributed Optimization of Self-concordant Empirical Loss, pages 289–341.  Springer International Publishing, Cham, 2018.
  * [43] X. Q. Zhize Li, Dmitry Kovalev and P. Richtárik.  Acceleration for compressed gradient descent in distributed and federated optimization.  International Conference on Machine Learning, 37, 2020.

##  Appendix A Proofs

###  A.1 Assumptions

######  Assumption 1.

The function $F$ is twice continuously differentiable.

######  Assumption 2.

There exist positive constants $\mu$ and $L$ such that $\mu
I\preceq\nabla^{2}F(w)\preceq LI,$ $\forall w\in\mathbb{R}^{d}.$

######  Assumption 3.

The function $F$ is bounded below by a scalar $\hat{F}$.

######  Assumption 4.

The gradients of $F$ are $L$-Lipschitz continuous for all
$w\in\mathbb{R}^{n}$.

###  A.2 Proof of Lemma 1

######  Lemma 1.

Let Assumption 1 hold. Let $A_{k}$ be defined as in (10) depending on the
choice of iterates update. Then, for any $k\geq 0$ there exist constants
$0<\mu_{1}\leq\mu_{2}$ such that

$\mu_{1}I\preceq A_{k}\preceq\mu_{2}I.$

###### Proof.

Truncated L-SR1 update. In this case
$A_{k}=(|B_{k+1}|_{\omega}^{\Omega})^{-1}$. By spectral decomposition
$B_{k}=V_{k}\Lambda_{k}V_{k}^{T}$. By truncation (Definition 1)

$\omega I\preceq|\Lambda_{k}|_{\omega}^{\Omega}\preceq\Omega I.$

Furthermore, $\dfrac{1}{\Omega}I\preceq A_{k}\preceq\dfrac{1}{\omega}I$.
Defining $\mu_{1}:=\frac{1}{\Omega}$ and $\mu_{2}:=\frac{1}{\omega}$, and
noting that $0<\mu_{1}\leq\mu_{2}$, gives the result.

FedSONIA update. In this case
$A_{k}=\left(|B_{k+1}^{\text{SONIA}}|_{\omega}^{\Omega}\right)^{-1}+\rho_{k}(I-\widetilde{V}_{k}\widetilde{V}_{k}^{T})$.
By (12)

$A_{k}=\widetilde{V}_{k}(|\Lambda_{k}|_{\omega}^{\Omega})^{-1}\widetilde{V}_{k}^{T}=\widetilde{V}_{k}(|\Lambda_{k}|_{\omega}^{\Omega})^{-1}-\rho_{k}I)\widetilde{V}_{k}^{T}+\rho_{k}I\stackrel{{\scriptstyle\eqref{eq:search_direction_sonia}}}{{\succeq}}\rho_{k}I.$

By truncation (Definition 1) $\omega
I\preceq|\Lambda_{k}|_{\omega}^{\Omega}\preceq\Omega I.$ Therefore,
$A_{k}\preceq\frac{1}{\omega}I+\rho_{k}I\preceq\frac{2}{\omega}$. Defining
$\mu_{1}:=\frac{1}{\Omega}$ and $\mu_{2}:=\frac{2}{\omega}$, and noting that
$0<\mu_{1}\leq\mu_{2}$, gives the result.

∎

###  A.3 Proof of Theorem 1

######  Theorem 1.

Let Assumptions 1, and 2 hold. Then for iterates $\\{w_{k}\\}$ generated by
Algorithm 5 with $0<\alpha_{k}=\alpha\leq\frac{\mu_{1}}{\mu_{2}^{2}L}$ we have

$F(w_{k})-F^{\star}\leq(1-\alpha\mu\mu_{1})^{k}[F(w_{0})-F^{\star}]$

for all $k\geq 0$.

###### Proof.

We have

$\displaystyle F(w_{k+1})=F(w_{k}-\alpha A_{k}\nabla F(w_{k}))$
$\displaystyle\stackrel{{\scriptstyle\text{Assump.
\ref{assum:strong_conv}}}}{{\leq}}F(w_{k})+\nabla F(w_{k})^{T}(-\alpha
A_{k}\nabla F(w_{k}))+\frac{L}{2}\|\alpha A_{k}\nabla F(w_{k})\|^{2}$
$\displaystyle\stackrel{{\scriptstyle
Lem.\ref{lem:positive_df}}}{{\leq}}F(w_{k})-\alpha\mu_{1}\|\nabla
F(w_{k})\|^{2}+\frac{\alpha^{2}\mu_{2}^{2}L}{2}\|\nabla F(w_{k})\|^{2}$
$\displaystyle=F(w_{k})-\alpha\left(\mu_{1}-\alpha\frac{\mu_{2}^{2}L}{2}\right)\|\nabla
F(w_{k})\|^{2}{\leq}F(w_{k})-\alpha\frac{\mu_{1}}{2}\|\nabla F(w_{k})\|^{2},$
(15)

where the last inequality is due to the choice of the step-size. By Assumtion
2, we have $2\mu(F(w)-F^{\star})\leq\|\nabla F(w)\|^{2}$. Therefore,

$\displaystyle F(w_{k+1})\leq F(w_{k})-\alpha\mu\mu_{1}(F(w_{k})-F^{\star}).$

Subtracting $F^{\star}$ from both sides,

$\displaystyle
F(w_{k+1})-F^{\star}\leq(1-\alpha\mu\mu_{1})(F(w_{k})-F^{\star}).$

Recursive application of the above inequality yields the desired result. ∎

###  A.4 Proof of Theorem 2

We begin with the technical lemma for FLECS with Truncated Hessian inverse
Approximation step (Algorithm 3).

######  Lemma 2.

Let Assumptions 1, 2 hold and $0\leq\omega\leq\mu$, $\Omega\geq L$. Then for
Hessian approximations $B_{k+1}$ we have

$\||B_{k+1}|_{\omega}^{\Omega}-\nabla^{2}f(w^{*})\|_{\text{F}}\leq\|B_{k+1}-\nabla^{2}f(w^{*})\|_{\text{F}}$

for all $k\geq 0$.

###### Proof.

By Definition 1

$|B_{k+1}|_{\omega}^{\Omega}\
:=V_{k+1}|\Lambda_{k}|_{\omega}^{\Omega}V_{k+1}^{T},$

where $B_{k+1}=V_{k+1}\Lambda_{k}V_{k+1}^{T}$ is spectral decomposition of
$B_{k+1}$. Then,

$\displaystyle\left\||B_{k+1}|_{\omega}^{\Omega}-\nabla^{2}f(w^{*})\right\|_{\text{F}}=\left\|V_{k+1}|\Lambda_{k+1}|_{\omega}^{\Omega}V_{k+1}^{T}-\nabla^{2}f(w^{*})\right\|_{\text{F}}$
$\displaystyle\stackrel{{\scriptstyle
V_{k+1}^{T}V_{k+1}=I}}{{=}}\left\|V_{k+1}|\Lambda_{k+1}|_{\omega}^{\Omega}V_{k+1}^{T}-V_{k+1}V_{k+1}^{T}\nabla^{2}f(w^{*})V_{k+1}V_{k+1}^{T}\right\|_{\text{F}}$
$\displaystyle=\left\||\Lambda_{k+1}|_{\omega}^{\Omega}-V_{k+1}^{T}\nabla^{2}f(w^{*})V_{k+1}\right\|_{\text{F}}.$

Analogically,

$\left\|B_{k+1}-\nabla^{2}f(w^{*})\right\|_{\text{F}}=\left\|\Lambda_{k}-V_{k+1}^{T}\nabla^{2}f(w^{*})V_{k+1}\right\|_{\text{F}}.$

Now, let us denote $\lambda_{i},~{}i=1,\ldots,d$ as eigenvalues of $B_{k+1}$
(diagonal elements of $\Lambda_{k+1}$). Let
$\widetilde{\lambda}_{1},\widetilde{\lambda}_{2},\ldots,\widetilde{\lambda}_{n}$
be corresponding diagonal elements of $|\Lambda_{k+1}|_{\omega}^{\Omega}$.

Noting, that

$0\preceq\omega I\preceq\mu I\preceq\nabla^{2}F(w^{*})\preceq LI\preceq\Omega
I,$

we have $|\widetilde{\lambda}_{i}-h_{i}|\leq|\lambda_{i}-h_{i}|$ for all
$i=1,\ldots,n$, where
$h_{i}=\left[V_{k+1}^{T}\nabla^{2}F(w^{*})V_{k+1}\right]_{ii}$.

Therefore,

$\||B_{k+1}|_{\omega}^{\Omega}-\nabla^{2}f(w^{*})\|_{\text{F}}\leq\|B_{k+1}-\nabla^{2}f(w^{*})\|_{\text{F}}$

by definition of Frobenius norm and the fact that matrices
$|B_{k+1}|_{\omega}^{\Omega}-\nabla^{2}f(w^{*})$ and
$\Lambda_{k}-V_{k+1}^{T}\nabla^{2}f(w^{*})V_{k+1}$ differ only in diagonal
elements. ∎

######  Theorem 2.

Let Assumptions 1 and 2 hold, and $\|w^{0}-w^{*}\|\leq\frac{\mu^{2}}{2}$,
$\frac{1}{n}\sum\limits_{i=1}^{n}\|B_{k+1}^{i}-\nabla^{2}f_{i}(w^{*})\|_{\text{F}}\leq\frac{2\mu^{2}}{L^{2}}$.
Then for iterates $\\{w_{k}\\}$ generated by Algorithm 5 with Truncated
inverse Hessian approximation step (Algorithm 3) with parameters
$\alpha_{k}=\alpha=1,~{}0\leq\omega\leq\mu,~{}\Omega\geq L$ we have

$\|w_{k}-w^{*}\|^{2}\leq\tfrac{1}{2^{k}}\|w_{0}-w^{*}\|^{2}$

for all $k\geq 0$.

###### Proof.

$\displaystyle\|w_{k+1}-w^{*}\|^{2}=\|w_{k}-w^{*}-A_{k}\nabla
F(w_{k})\|^{2}\leq\|A_{k}\|^{2}\|A_{k}^{-1}(w_{k}-w^{*})-\nabla
F(w_{k})\|^{2}$ $\displaystyle\stackrel{{\scriptstyle\text{Lem.
\ref{lem:positive_df}}}}{{\leq}}\mu_{2}^{2}\|A_{k}^{-1}(w_{k}-w^{*})-\nabla
F(w_{k})\|^{2}$ $\displaystyle\leq
2\mu_{2}^{2}\left(\|\left(A_{k}^{-1}-\nabla^{2}F(w^{*})\right)(w_{k}-w^{*})\|^{2}+\|\nabla^{2}F(w^{*})(w_{k}-w^{*})-\nabla
F(w_{k})+\nabla F(w^{*})\|^{2}\right)$
$\displaystyle\stackrel{{\scriptstyle\text{Assump.
\ref{assum:strong_conv}}}}{{\leq}}2\mu_{2}^{2}\left(\|A_{k}^{-1}-\nabla^{2}F(w^{*})\|^{2}\|w_{k}-w^{*}\|^{2}+\frac{L^{2}}{4}\|w_{k}-w^{*}\|^{4}\right)$
$\displaystyle\leq
2\mu_{2}^{2}\|w_{k}-w^{*}\|^{2}\left(\|A_{k}^{-1}-\nabla^{2}F(w^{*})\|_{\text{F}}^{2}+\frac{L^{2}}{4}\|w_{k}-w^{*}\|_{\text{F}}^{2}\right)$
$\displaystyle\stackrel{{\scriptstyle\eqref{eq:search_dir_sr1}}}{{=}}\leq
2\mu_{2}^{2}\|w_{k}-w^{*}\|^{2}\left(\||B_{k+1}|_{\omega}^{\Omega}-\nabla^{2}F(w^{*})\|_{\text{F}}^{2}+\frac{L^{2}}{4}\|w_{k}-w^{*}\|_{\text{F}}^{2}\right)$
$\displaystyle\stackrel{{\scriptstyle\text{Lem.\ref{lem:trunc}}}}{{\leq}}2\mu_{2}^{2}\|w_{k}-w^{*}\|^{2}\left(\|B_{k+1}-\nabla^{2}F(w^{*})\|_{F}^{2}+\frac{L^{2}}{4}\|w_{k}-w^{*}\|^{2}\right)$

Therefore, by assumptions of the theorem

$\displaystyle\|w_{k+1}-w^{*}\|^{2}\leq\frac{1}{\mu^{2}}\|w_{k}-w^{*}\|^{2}\left(\|B_{k+1}-\nabla^{2}F(w^{*})\|_{\text{F}}^{2}+\frac{L^{2}}{4}\|w_{k}-w^{*}\|^{2}\right)\leq\frac{1}{2}\|w_{k}-w^{*}\|^{2}$

∎

###  A.5 Proof of Theorem 3

######  Theorem 3.

Let Assumptions1, 3, 4 hold and $w_{0}$ be the starting point. Then after
$T>0$ iterations of Algorithm 5 with
$0<\alpha_{k}=\alpha\leq\frac{\mu_{1}}{\mu_{2}^{2}L}$, we have

$\tfrac{1}{T}\textstyle{\sum}_{k=0}^{T-1}\|\nabla
F(w_{k})\|^{2}\leq\tfrac{2[F(w_{0})-\hat{F}]}{\alpha\mu_{1}T}\stackrel{{\scriptstyle
T\to\infty}}{{\longrightarrow}}0.$

###### Proof.

By (15)

$\displaystyle F(w_{k+1})\leq F(w_{k})-\alpha\frac{\mu_{1}}{2}\|\nabla
F(w_{k})\|^{2}.$

Summing both sides of the above inequality from $k=0$ to $T-1$, we obtain

$\displaystyle\widehat{F}-F(w_{0})\stackrel{{\scriptstyle\text{Assum.
\ref{assum:boundness}}}}{{\leq}}F(w_{T})-F(w_{0})=\sum_{k=0}^{T-1}(F(w_{k+1})-F(w_{k}))\leq-\sum_{k=0}^{T-1}\alpha\frac{\mu_{1}}{2}\|\nabla
F(w_{k})\|^{2}.$

Therefore, we have

$\displaystyle\sum_{k=0}^{T-1}\|\nabla F(w_{k})\|^{2}$
$\displaystyle\leq\frac{2[F(w_{0})-\widehat{F}]}{\alpha\mu_{1}}.$ (16)

Dividing (16) by $T$ we conclude

$\displaystyle\frac{1}{T}\sum_{k=0}^{T-1}\|\nabla F(w_{k})\|^{2}$
$\displaystyle\leq\frac{2[F(w_{0})-\widehat{F}]}{\alpha\mu_{1}T}.$

∎

##  Appendix B Additional Numerical Experiments

###  B.1 Top-$K$ Compressor

In this experiment we used Top-$K$ compression with $k=4d$ on a9a dataset with
regularization parameter $\mu=10^{-5}$. We compare FedNL [36] and FLECS in the
number of hessian-vector products, since compression cost is roughly the same.
For FLECS we set memory size $m=16$ and step-size $\alpha=1$. For both
algorithms we set $B_{0}^{i}=\nabla F(w_{0})$ and Hessian approximation rate
$\beta=1$. In this experiment we omit initialization cost (computation of the
full Hessian at starting point) for both algorithms.

Figure 2: Comparison of objective function $F(w_{k})$ and the squared norm of
gradient $\|\nabla F(w_{k})\|^{2}$ between FLECS and FedNL.

###  B.2 Comparison with DIANA and ADIANA.

In order to compare our framework with DIANA and ADIANA we use random
dithering with number of levels $s=\sqrt{d}$ and $p=\infty$. For FLECS we set
hyperparametes as: $m=1$, $\beta=1$, $\alpha=1$ and initialize Hessian
approximations with $0$.

For DIANA and ADIANA we used theoretical parameters $L=1/4$ and strong
convexity constant $\mu=10^{-5}$. For these methods we used an original code
provided by authors [27, 43].

Figure 3: Comparison of objective function $F(w_{k})$ and the squared norm of
gradient $\|\nabla F(w_{k})\|^{2}$ between DIANA, ADIANA and FLECS.

###  B.3 Top-$K$ vs Random Dithering

In this subsection we want to compare Top-$K$ compressor with random dithering
for FLECS. We set hyperparametes as: $m=16$, $\beta=1$, $\alpha=1$ and
initialize Hessian approximations with $\nabla^{2}f(w_{0})$. We use Top-$K$
with $k=4d$ and random dithering with $s=128$ and $p=\infty$, so communication
cost per iteration is roughly the same. In this experiment we omit
initialization cost (computation of the full Hessian at starting point) for
both algorithms.

Figure 4: Comparison of objective function $F(w_{k})$ and the squared norm of
gradient $\|\nabla F(w_{k})\|^{2}$ .

###  B.4 Dependence on the Memory Size

Here we compare different memory sizes for FLECS on a9a dataset. We use random
dithering with parameters $s=128$ and $p=\infty$. Hyperparametes set as
$\beta=1$, $\alpha=1$ and initialize Hessian approximations initialized with
$0$.

Figure 5: Comparison of objective function $F(w_{k})$ and the squared norm of
gradient $\|\nabla F(w_{k})\|^{2}$ for different memory sizes.

##  Appendix C Compressors

FLECS allows several different compression operators. In this section we give
definitions and provide examples of possible compressors for our framework.

###  C.1 Unbiased Compressors.

######  Definition 2.

A random operator $\mathcal{C}:~{}\mathbb{R}^{d\times m}\to\mathbb{R}^{d\times
m}$ satisfying

${\mathbb{E}}_{\mathcal{Q}}\left[\mathcal{C}(X)\right]=X,~{}~{}{\mathbb{E}}_{\mathcal{C}}\left[\|\mathcal{C}(X)\|_{\text{F}}^{2}\right]\leq(\omega+1)\|X\|_{\text{F}}^{2}],$
(17)

for all $X\in\mathbb{R}^{d\times m}$ is a unbiased compressor operator.

Random dithering for vectors ($\omega$-quantization) [2], applied to each
column of the matrix, is an example of unbiased compressors that we used in
our experiments.

###  C.2 Biased compressors.

######  Definition 3.

A random operator $\mathcal{C}:~{}\mathbb{R}^{d\times m}\to\mathbb{R}^{d\times
m}$ satisfying

${\mathbb{E}}_{\mathcal{C}}\left[\mathcal{C}(X)\right]\leq\|X\|_{\text{F}},~{}~{}{\mathbb{E}}_{\mathcal{C}}\left[\|\mathcal{C}(X)\|_{\text{F}}^{2}\right]\leq(1-\delta)\|X\|_{\text{F}}^{2}],$
(18)

for all $X\in\mathbb{R}^{d\times m}$ is a contractive compressor operator.

In our experiments we use Top-$K$ as a contractive compressor.