# A test of Radial Acceleration Relation for the Giles et al Chandra cluster
sample

S. Pradyumna<EMAIL_ADDRESS>Shantanu Desai  E-mail:
<EMAIL_ADDRESS>Department of Physics, Indian Institute of Technology,
Hyderabad, Telangana-502285, India

###### Abstract

We carry out a test of the radial acceleration relation (RAR) for a sample of
10 dynamically relaxed and cool-core galaxy clusters imaged by the Chandra
X-ray telescope, which was studied in Giles et al. (2017). For this sample, we
observe that the best-fit RAR shows a very tight residual scatter equal to
0.09 dex. We obtain an acceleration scale of $1.59\times 10^{-9}m/s^{2}$,
which is about an order of magnitude higher than that obtained for galaxies.
Furthermore, the best-fit RAR parameters differ from those estimated from some
of the previously analyzed cluster samples, which indicates that the
acceleration scale found from the RAR could be of an emergent nature, instead
of a fundamental universal scale.

##  I Introduction

Three recent works Tian et al. (2020); Chan and Del Popolo (2020); Pradyumna
et al. (2021) have studied disparate samples of galaxy clusters, in order to
ascertain if they obey the Radial Acceleration relation (RAR). The RAR is a
tight correlation between the baryonic ($a_{bar}$) and the total dynamical
acceleration ($a_{tot}$). A first definitive case for the RAR was asserted
using the SPARC sample consisting of spiral galaxies McGaugh et al. (2016),
who showed the RAR is obeyed with a scatter of 0.13 dex, with most of this
scatter been attributed to observational uncertainties Li et al. (2018). The
RAR, which was deduced from the SPARC sample can be written as follows:

$a_{tot}=\frac{a_{bar}}{1-e^{-\sqrt{a_{bar}/a_{0}}}},$ (1)

where $a_{0}\sim 1.2\times 10^{-10}m/s^{2}$ McGaugh et al. (2016). However,
some objections to $a_{0}$ been a fundamental universal constant (valid for
all galaxies) have been raised Rodrigues et al. (2018a, b); Rodrigues et al.
(2020); Marra et al. (2020). Ref. Rodrigues et al. (2018a) showed that the
presence of a fundamental acceleration scale is ruled out at more than
10$\sigma$, and there is no value of $a_{0}$ which is compatible within
5$\sigma$ for all galaxies. This work concluded that the RAR is an emergent
one, which arises after stacking the data from a large number of galaxies.
There is an ongoing debate on these issues McGaugh et al. (2018); Kroupa et
al. (2018); Rodrigues et al. (2018b); Rodrigues et al. (2020); Marra et al.
(2020). Most recently, the results in Rodrigues et al. (2018a) were also
confirmed in Chang and Zhou (2019); Marra et al. (2020). (See also Ref. Zhou
et al. (2020)).

More generally, the RAR can be re-written as a linear regression relation
between the logarithms of $a_{bar}$ and $a_{tot}$,
$\ln(a_{tot})=m\ln(a_{bar})+b$. All RAR studies for clusters have been
obtained by carrying out this linear regression between these two
accelerations in logarithmic space Tian et al. (2020). This relation can be
trivially deduced from Milgrom’s MOND paradigm Famaey and McGaugh (2012).
However, the observed tight scatter cannot be trivially predicted by the
standard $\Lambda$CDM model Desmond (2017); Navarro et al. (2017); Ludlow et
al. (2017) (see also the discussion in McGaugh (2020)). However, if $a_{0}$ is
of emergent nature, as argued in Rodrigues et al. (2018a); Marra et al.
(2020), then the RAR could be reproduced from the complex interplay between
dark matter and baryons Ludlow et al. (2017); Navarro et al. (2017); Dutton et
al. (2019), including the effects due to the adiabatic contraction of dark in
the presence of baryons Paranjape and Sheth (2021).

We have been knowing for more than three decades that MOND does not work for
relaxed galaxy clusters The and White (1988); Seeram and Desai (2021) (and
references therein). Furthermore, relativistic MOND theories cannot explain
the coincident GW-EM signal Boran et al. (2018) as well as large-scale
structure probes Pardo and Spergel (2020). However, a detailed
characterization of the RAR in galaxy clusters could be an invaluable probe,
that could be used to adjudicate between various alternatives to the standard
$\Lambda$CDM model, which can reproduce the success of Milgrom’s MOND’s
paradigm, but are consistent with $\Lambda$CDM predictions at larger scales
from CMB Ade et al. (2016).

Motivated by these considerations, three different groups recently carried out
a test of the RAR for multiple galaxy cluster samples Tian et al. (2020); Chan
and Del Popolo (2020); Pradyumna et al. (2021). They have found that the
residual scatter between the observed data and best-fit RAR is between
0.11-0.13 dex (see Table I of Pradyumna et al. (2021) for a summary). However,
the aforementioned works have found that the acceleration scale $a_{0}$ in Eq.
1 is $\mathcal{O}(10^{-9})m/s^{2}$, which is about an order of magnitude
larger than that obtained for rotationally supported galaxies McGaugh et al.
(2016). Most recently, this work has also been extended to galaxy groups
Gopika and Desai (2021), who find that the scatter is larger, and $a_{0}$
falls in-between single galaxies and clusters.

Given the importance of RAR as a probe of $\Lambda$CDM and its alternatives,
it is important to cross-check these results with more cluster samples to gain
further insight on how universal this relation is for clusters, and whether
there are any exceptions. In the past three decades a large number of
dedicated X-ray based cluster surveys using ROSAT, Chandra, and XMM-Newton
(along with other multi-wavelength data) have used the observed cluster
samples to address a plethora of science goals in Cosmology, galaxy evolution,
and fundamental Physics Allen et al. (2011); Vikhlinin et al. (2014). The same
data can also be used for precision tests of RAR.

In this work, we therefore carry out another test of RAR with one such galaxy
cluster sample studied in Giles et al. (2017), which we refer to as Giles et
al cluster sample. This sample has been recently used to test Verlinde’s
Emergent gravity paradigm Halenka and Miller (2020). The outline of this
manuscript is as follows. The basic methodology used for our analysis is
described in Sect. II. The Giles et al cluster sample along with our data
analysis procedure is outlined in Sect. III. Our results are described in
Sect. IV. We conclude in Sect. V.

##  II Basic Cluster Physics

We now provided an abridged description of all the ingredients necessary to
carry out a test of RAR using clusters. More details can also be found in our
past works Gopika and Desai (2020); Pradyumna et al. (2021); Gopika and Desai
(2021).

Galaxy clusters are the largest gravitationally collapsed objects in the
universe Vikhlinin et al. (2014). Most of their mass ($\sim 85\%$) comprises
dark matter, whereas most of the remaining baryonic mass resides in the form
of diffuse hot gas ($\sim 13\%$) in the Intra-Cluster Medium (ICM). Only a
small fraction of the baryonic mass ($\sim 2\%$) is in the form of stars. At
X-Ray wavelengths, the hot ICM emits via thermal bremsstrahlung Allen et al.
(2011). The total mass of a cluster can be estimated from X-Ray measurements
assuming that the ICM is in hydrostatic equilibrium within the gravitational
potential of the cluster. The total acceleration for test particles in
hydrostatic equilibrium within a gravitational potential is given by Sarazin
(1986); Allen et al. (2011):

$a_{tot}(r)=-\frac{k_{b}T}{r\mu m_{p}}\left(\frac{d\ln\rho_{gas}}{d\ln
r}+\frac{d\ln T}{d\ln r}\right),$ (2)

Here, $a_{tot}$ is the dynamic acceleration of the cluster as a function of
radius from the cluster center ($r$), $k_{b}$ is the Boltzmann constant, $T$
is the parametrized temperature, $\rho_{gas}$ is the parametrized density of
the gas, $\mu m_{p}$ is the mean molecular weight. We use a $\mu$ value of 0.6
resulting in $\mu m_{p}\approx 10^{-27}\,\mathrm{kg}$. The hydrostatic
equilibrium equation is applicable for relaxed clusters. Once we know the
total mass we can estimate the total acceleration $(a_{t})$ from Newtonian
gravity.

Now, for estimating the acceleration due to baryons, we need to estimate the
gas mass ($M_{gas}$) and the stellar mass ($M_{star}$). The gas density
($\rho_{gas}$) profiles are estimated and modelled using the observed
projected X-ray emissivity profile. Assuming spherical symmetry, we can derive
the gas mass $M_{gas}$ using:

$M_{gas}=\int 4\pi r^{2}\rho_{gas}(r)dr$ (3)

It should be noted that the aforementioned calculations of gas and total mass
assume spherical symmetry. The uncertainty associated with this assumption is
about 5% (See Gopika and Desai (2020) and references therein). For stellar
mass, we use the empirical relation obtained in Lin et al. (2012), which gives
the stellar mass at $R_{500}$, where $R_{500}$ is the radius at which the
overdensity is 500 times the critical density of the universe, i.e.
$\rho(R_{500})=500\left(\frac{3H^{2}}{8\pi G}\right)$. The total star mass
$M_{star}(R_{500})$ at $R_{500}$ is given by:

$\frac{M_{star}(R_{500})}{10^{12}M_{\odot}}=1.8\left(\frac{M_{500}}{10^{14}M_{\odot}}\right)^{0.71}$
(4)

Then, to estimate the stellar mass at any other radius ($r$), we assume that
the stellar component follows an isothermal distribution, and we thereby
obtain Rahvar and Mashhoon (2014):

$M_{star}(r)=\left(\frac{r}{R_{500}}\right)M_{star}(R_{500})$ (5)

The baryonic acceleration ($a_{b}$) can then be obtained at any value of the
radius $r$ from the sum of the gas and star mass.

##  III Data and Analysis:

Giles et al Giles et al. (2017) (G17, hereafter) provide the temperature and
gas density profile fits for 34 clusters within the redshift range $0.15\leq
z\leq 0.3$, which were imaged with the Chandra X-ray telescope. The main goal
of this work was to characterize the scatter in the Luminosity-Mass relation
for this sample. This cluster sample was selected from the ROSAT based
Brightest cluster sample survey Ebeling et al. (1998, 2000) with a lower
luminosity cutoff of $L_{X}$ (0.1-2.4 keV) $\geq 6\times 10^{44}$ erg/sec
Dahle (2006). This cut provided a sample of 36 clusters. From these, A689 was
culled because of data quality issues and Zw5768 did not pass the luminosity
cut using a revised redshift estimate. The exposure times for these clusters
are in the range of 7.3-450 ks (cf. Table 1 of G17). The luminosity of these
clusters in X-Ray band of 0.1-2.4 keV falls between $\sim(6.3-24.5)\times
10^{44}$ erg/s. More details of the observations and data reduction can be
found in G17.

Since hydrostatic equilibrium can be robustly used to obtain the mass only if
the cluster is dynamically relaxed and has a cool core, we select a sub-sample
of relaxed as well as cool core clusters for our analysis. This was done in
two steps: (1) the clusters were determined to be in a relaxed or non-relaxed
state, on the basis of the centroid shift ($\langle w\rangle$) Poole et al.
(2006), the description of which can be found G17; (2) the cluster sample was
further bifurcated into cool core and non-cool core sample. A cut of $\langle
w\rangle<0.009r_{500}$ was made to select dynamically relaxed clusters. The
presence of a cool core was based on the calculating the central cooling time
given by:

$t_{cool}(yr)=8.5\times
10^{10}(\frac{n_{p}}{10^{-3}})^{-1}(0.079kT_{CCT})^{1/2}$

The last parameter used for the selection of cool core sample was the
cuspiness, which is defined as the logarithmic slope of gas density at
$0.04r_{500}$ Hudson et al. (2010). In summary, cuts were applied on $\langle
w\rangle<0.009r_{500}$, $t_{cool}<7.7Gyr$, cuspiness $>0.7$, to select the
dynamically relaxed and cool core sample. After applying these cuts, we are
left with a sample of 10 clusters (from the original dataset of 34 clusters)
for our analysis. They are: A2204, RXJ1720.1+2638, A1423, Z2089,
RXJ2129.6+0005, A1835, MS1455.0+2232, RXJ0439.0+0715, RXJ0437.1+0043, Z3146.

The gas density ($\rho_{gas}$) model used is equal to
$1.624m_{p}\sqrt{n_{p}(r)n_{e}(r)}$, where $m_{p}$ is the proton mass; $n_{p}$
and $n_{e}$ are the proton and electron number densities, respectively. Their
product ($n_{e}n_{p}$) is parametrized by a modified single-$\beta$ profile:

$n_{e}(r)n_{p}(r)=\frac{(r/r_{c})^{-\alpha}}{(1+r^{2}/r_{c}^{2})^{3\beta-\alpha/2}}\frac{n_{0}^{2}}{(1+r^{\gamma}/r_{s}^{\gamma})^{\epsilon/\gamma}}$
(6)

This is a variant of the double-$\beta$ profile used by Vikhlinin et al
Vikhlinin et al. (2006), which was also used in our previous works Gupta and
Desai (2019, 2020); Gopika and Desai (2020); Pradyumna et al. (2021) The
second term in the double-$\beta$ profile was dropped so as to fit both the
low and high quality data in the sample.

The temperature profile used for our current analysis is:

$T(r)=T_{0}\frac{(x+T_{min}/T_{0})}{x+1}\frac{(r/r_{t})^{-a}}{\left[1+(r/r_{t})^{b}\right]^{c/b}},$
(7)

where $x=\left(\frac{r}{r_{cool}}\right)^{a_{cool}}$. The values of the free
parameters in Eq. 7 for these clusters can be found in Table A2 of G17. The
error estimates for the temperature and density have been obtained in G17
using Monte-Carlo simulations. For this work, the temperature error bars used
in our analysis have been obtained by digitization of the plots in G17, which
in turn were used to propagate the errors in $M_{tot}$. For the errors on the
gas density, we again digitized the emission measure plots given in G17. We
calculated the fractional errors in the emission measure data at the radii of
interest by using spline interpolation. Then using the fact that the gas
density is proportional to the square root of the emission measure, the
fractional error in the density is given as:
$\frac{\Delta\rho_{gas}}{\rho_{gas}}=\frac{1}{2}\frac{\Delta EM}{EM}$ where EM
stands for the emission measure value. We now plug in the values for $\rho(r)$
and $T(r)$ from Eqns. 6 and 7 respectively in Eq. 2 to calculate the total
dynamical acceleration $a_{tot}(r)$.

To estimate baryonic mass, we need to determine the gas and star mass. The gas
mass can be estimated from the gas density profile in Eq. 6. To estimate the
stellar mass, we need $R_{500}$ and $M_{500}$ for these clusters. G17 provides
these $R_{500}$ estimated using the hydrostatic mass estimates. We plugged
these values in Eq. 4 to obtain the total stellar mass estimate at $R_{500}$.
The stellar mass at all other radii was obtained from Eq. 5. Therefore, the
baryonic acceleration can then estimated from this total baryonic mass
assuming Newtonian gravity.

Now, we do a joint fit to the RAR using the calculated data along with
associated errors for $a_{tot}$ and $a_{bar}$ at $r$ = 100, 200, 400, and 600
kpc for all the clusters. The only exception is Z2089, for which we used
measurements at 100, 200, 400, and 560 kpc. Similar to previous works Tian et
al. (2020); Pradyumna et al. (2021), we do a fit in log-acceleration space by
doing a linear regression to the relation:

$y=mx+b$ (8)

where $y=\ln(a_{tot})$ and $x=\ln(a_{bar})$. We then maximize a log-likelihood
given by :

$\displaystyle-2\ln L$ $\displaystyle=$ $\displaystyle\large{\sum_{i}}\ln
2\pi\sigma_{i}^{2}+\large{\sum_{i}}\frac{[y_{i}-(mx_{i}+b)]^{2}}{\sigma_{i}^{2}}$
(9) $\displaystyle\sigma_{i}^{2}$ $\displaystyle=$
$\displaystyle\sigma_{y_{i}}^{2}+m^{2}\sigma_{x_{i}}^{2}+\sigma_{int}^{2}$
(10)

Here, $\sigma_{y_{i}}$ and $\sigma_{x_{i}}$ denote the errors in $y$ and $x$
respectively; $\sigma_{int}$ is the intrinsic scatter for our linear relation
between the logarithms of the accelerations, which is a free parameter. We
used the emcee MCMC sampler Foreman-Mackey et al. (2013) to sample the
likelihood in Eq. 9. The total number of steps and burn-in steps used for the
MCMC analysis are 3000 and 500 respectively. The number of walkers used for
our MCMC analysis is equal to 40. We used uniform non-informative priors on
$m$, $b$, and $\ln(\sigma_{int})$ given by: $m\in[-10,10]$ ; $b\in[-100,100]$;
and $\ln(\sigma_{int})\in[-5,-1]$.

For calculating the residual scatter in dex, we first calculate residuals by
subtracting the $\log_{10}(a_{dyn})$ predicted by our best-fit from
$\log_{10}(a_{dyn})$ calculated from the data. This is the same practise used
in the literature to characterize the scatter in the RAR Pradyumna et al.
(2021). We then fit a normal distribution with a mean of zero to these
residuals. The standard deviation obtained is used to characterize the
residual scatter in dex.

##  IV Results

We now present the results of our analysis of 10 relaxed cool core clusters
from this sample. Fig. 1 shows the best-fit line obtained for the data along
with the best fits obtained by McGaugh et al McGaugh et al. (2016) (for
galaxies), Tian et al Tian et al. (2020), Chan and Del Popolo Chan and Del
Popolo (2020) and our previous analysis Pradyumna et al. (2021) (on galaxy
clusters). Fig 2 shows the 68% and 95% credible intervals for the slope,
intercept, and the natural log of intrinsic scatter obtained using the
ChainConsumer post-processing package Hinton et al. (2020), using the MCMC
chains from emcee as inputs. The best-fit slope and intercept obtained are
$m=1.03\pm 0.05$ and $b=2.7^{+1.3}_{-1.2}$ respectively. The residual scatter
obtained for the best fit is about 0.08 dex. The best-fit values and the
scatter are shown in Table 1 alongside the fits obtained in previous works for
comparison. We note that the best-fit slope and intercept observed for XCOP
sample and the current sample agree within $1\sigma$. However, the best-fit
parameters of the G17 sample differ from the other cluster samples, such as
the Vikhlinin et al Chandra sample Vikhlinin et al. (2006), HIFLUGS sample
Chan and Del Popolo (2020), and the CLASH sample Tian et al. (2020) to between
$3-5\sigma$. The acceleration scale that we obtain for current sample on
fixing the slope to $m=0.5$ is equal to $(1.59\pm 0.14)\times 10^{-9}m/s^{2}$.
The error in $a_{0}$ was estimated using error propagation, based on the
$1\sigma$ error in $b$, which was obtained from redoing the RAR fit for a
fixed slope of 0.5. The acceleration scale obtained for this sample of galaxy
clusters is therefore, an order of magnitude higher than the $a_{0}$ obtained
from the SPARC sample McGaugh et al. (2016) for galaxies, in agreement with
all previous analyses carried out on galaxy clusters. Tian et al. (2020); Chan
and Del Popolo (2020); Pradyumna et al. (2021).

When we try to apply a cut to look for dynamically relaxed but non-cool
clusters, using only the dynamically relaxed criteria mentioned in Sec III, we
find three additional clusters for our analysis. Then, we tried to do a
similar analysis on all the relaxed clusters (both cool core and non-cool
core). The best fit slope and intercept obtained in this case are $m=0.92\pm
0.07$ and $b=0.0\pm 1.7$, and the scatter increases slightly to 0.095 dex.

Therefore, our results reinforce those in our previous work, and indicate the
RAR is obeyed for clusters with a tight intrinsic scatter. However, this
relation is not universal, as all the cluster samples do not agree on the RAR
best-fit parameters. These inconsistent fits among the different cluster
samples, reaffirm the view that $a_{0}$ could be of emergent nature, instead
of a fundamental acceleration scale, as previously discussed Rodrigues et al.
(2018a); Marra et al. (2020).

Cluster Sample | Slope | Intercept | Intrinsic scatter |  $a_{0}$ $(m/s^{2})$ | Scatter (dex)
---|---|---|---|---|---
Giles Chandra (this work) | $1.03^{+0.05}_{-0.05}$ | $2.7^{+1.3}_{-1.2}$ | $0.07_{-0.05}^{+0.03}$ | $(1.59\pm 0.14)\times 10^{-9}$ | 0.08
XCOP Sample Pradyumna et al. (2021) | $1.09^{+0.07}_{-0.07}$ | $4.44^{+1.80}_{-1.85}$ | $0.21_{-0.02}^{+0.03}$ | $(1.12\pm 0.11)\times 10^{-9}$ | 0.11
Vikhlinin Chandra Pradyumna et al. (2021) | $0.77\pm 0.1$ | $-3.5^{+2.6}_{-2.7}$ | ($0.0002\pm 0.018$)% | $(9.26\pm 1.66)\times 10^{-10}$ | 0.14
Chan et al Chan and Del Popolo (2020) | $0.72_{-0.06}^{+0.06}$ | $-5.49_{-1.39}^{+1.39}$ | $0.27_{-0.02}^{+0.02}$ | $9.5\times 10^{-10}$ | 0.18
Tian et al Tian et al. (2020) | $0.51^{+0.04}_{-0.05}$ | $-9.80^{+1.07}_{-1.08}$ | $14.7^{+2.9}_{-2.8}\%$ | $(2.02\pm 0.11)\times 10^{-9}$ | 0.11
Table 1: The best-fit values obtained for the Chandra sample of G17 along with
those obtained for other cluster samples. We see that the acceleration scale
we obtain from this sample is in agreement with results from previous RAR
analysis with clusters and about an order of magnitude higher than that
obtained for galaxies. Figure 1: The best-fit RAR(blue solid line) obtained
for $a_{tot}$ vs $a_{bar}$ data for a sample of 10 dynamically relaxed and
cool-core galaxy clusters from the G17 catalog. Best fits from previous works
on RAR for clusters have been overlaid for comparison. Figure 2: 68% and 95%
marginalized credible intervals obtained for the slope, intercept, and natural
log of the intrinsic scatter on fitting the RAR to the Giles et al Chandra
galaxy cluster sample Giles et al. (2017).

##  V Conclusions

We implement a test of the scatter in the correlation between the total
dynamical and baryonic acceleration in log-log space (dubbed as RAR) using a
sample of 10 dynamically relaxed cool core clusters, obtained from a dataset
of 35 ROSAT selected clusters imaged with Chandra, and described in detail in
G17. The aim of this work was to further assess the conclusions of previous
studies of RAR using ancillary cluster samples Tian et al. (2020); Pradyumna
et al. (2021); Chan and Del Popolo (2020). These previous analyses showed that
the intrinsic scatter of the RAR is the same as that seen for spiral galaxies
McGaugh et al. (2016), but the observed acceleration scale is elevated by an
order of magnitude, and the best-fit values are inconsistent among the cluster
samples Pradyumna et al. (2021).

The best-fit values for the parameters of the RAR using the G17 cluster sample
can be found in Fig. 2. The best-fit RAR along with the data is shown in Fig.
1. A comparison of our best-fit results along with those from the previous
works is summarized in Table 1. We find that the scatter for this sample is
equal to 0.08 dex and is slightly increased to 0.094 dex, once we include the
relaxed non-cool core clusters. Our estimate for the acceleration scale we
obtained by positing a slope of 0.5 in the RAR regression relation is about
about an order of magnitude larger than that deduced using the SPARC sample.
The best-fit value for the slope and intercept is in agreement with the same
for the XCOP sample, but differ from other cluster samples previously analyzed
by upto $5\sigma$. This shows that the RAR is not universal for clusters.
Therefore, this reinforces the fact alluded to in Rodrigues et al. (2018a);
Marra et al. (2020), that the acceleration scale deduced from RAR is not a
fundamental universal constant, but instead could be an emergent effect caused
by stacking data from different clusters.

###### Acknowledgements.

We would like to thank Paul Giles for useful correspondence about G17. We are
also grateful to the anonymous referee for useful feedback on the manuscript.

## References

  * Giles et al. (2017) P. A. Giles, B. J. Maughan, H. Dahle, M. Bonamente, D. Landry, C. Jones, M. Joy, S. S. Murray, and N. van der Pyl, Mon. Not. R. Astron. Soc. 465, 858 (2017), eprint 1510.04270.
  * Tian et al. (2020) Y. Tian, K. Umetsu, C.-M. Ko, M. Donahue, and I. N. Chiu, Astrophys. J. 896, 70 (2020), eprint 2001.08340.
  * Chan and Del Popolo (2020) M. H. Chan and A. Del Popolo, Mon. Not. R. Astron. Soc. p. 218 (2020), eprint 2001.06141.
  * Pradyumna et al. (2021) S. Pradyumna, S. Gupta, S. Seeram, and S. Desai, Physics of the Dark Universe 31, 100765 (2021), eprint 2011.06421.
  * McGaugh et al. (2016) S. S. McGaugh, F. Lelli, and J. M. Schombert, Phys. Rev. Lett. 117, 201101 (2016), eprint 1609.05917.
  * Li et al. (2018) P. Li, F. Lelli, S. McGaugh, and J. Schombert, Astron. & Astrophys. 615, A3 (2018), eprint 1803.00022.
  * Rodrigues et al. (2018a) D. C. Rodrigues, V. Marra, A. del Popolo, and Z. Davari, Nature Astronomy 2, 668 (2018a), eprint 1806.06803.
  * Rodrigues et al. (2018b) D. C. Rodrigues, V. Marra, A. Del Popolo, and Z. Davari, Nature Astronomy 2, 927 (2018b), eprint 1811.05882.
  * Rodrigues et al. (2020) D. C. Rodrigues, V. Marra, A. Del Popolo, and Z. Davari, Nature Astronomy 4, 134 (2020), eprint 2002.01970.
  * Marra et al. (2020) V. Marra, D. C. Rodrigues, and Á. O. F. de Almeida, Mon. Not. R. Astron. Soc. 494, 2875 (2020), eprint 2002.03946.
  * McGaugh et al. (2018) S. S. McGaugh, P. Li, F. Lelli, and J. M. Schombert, Nature Astronomy 2, 924 (2018).
  * Kroupa et al. (2018) P. Kroupa, I. Banik, H. Haghi, A. H. Zonoozi, J. Dabringhausen, B. Javanmardi, O. Müller, X. Wu, and H. Zhao, Nature Astronomy 2, 925 (2018), eprint 1811.11754.
  * Chang and Zhou (2019) Z. Chang and Y. Zhou, Mon. Not. R. Astron. Soc. 486, 1658 (2019), eprint 1812.05002.
  * Zhou et al. (2020) Y. Zhou, A. Del Popolo, and Z. Chang, Physics of the Dark Universe 28, 100468 (2020), eprint 2008.04065.
  * Famaey and McGaugh (2012) B. Famaey and S. S. McGaugh, Living Reviews in Relativity 15, 10 (2012), eprint 1112.3960.
  * Desmond (2017) H. Desmond, Mon. Not. R. Astron. Soc. 464, 4160 (2017), eprint 1607.01800.
  * Navarro et al. (2017) J. F. Navarro, A. Benítez-Llambay, A. Fattahi, C. S. Frenk, A. D. Ludlow, K. A. Oman, M. Schaller, and T. Theuns, Mon. Not. R. Astron. Soc. 471, 1841 (2017), eprint 1612.06329.
  * Ludlow et al. (2017) A. D. Ludlow, A. Benítez-Llambay, M. Schaller, T. Theuns, C. S. Frenk, R. Bower, J. Schaye, R. A. Crain, J. F. Navarro, A. Fattahi, et al., Phys. Rev. Lett. 118, 161103 (2017), eprint 1610.07663.
  * McGaugh (2020) S. McGaugh, Galaxies 8, 35 (2020), eprint 2004.14402.
  * Dutton et al. (2019) A. A. Dutton, A. V. Macciò, A. Obreja, and T. Buck, Mon. Not. R. Astron. Soc. 485, 1886 (2019), eprint 1902.06751.
  * Paranjape and Sheth (2021) A. Paranjape and R. K. Sheth, arXiv e-prints arXiv:2102.13116 (2021), eprint 2102.13116.
  * The and White (1988) L. S. The and S. D. M. White, Astron. J. 95, 1642 (1988).
  * Seeram and Desai (2021) S. Seeram and S. Desai, Journal of Astrophysics and Astronomy 42, 3 (2021).
  * Boran et al. (2018) S. Boran, S. Desai, E. O. Kahya, and R. P. Woodard, Phys. Rev. D 97, 041501 (2018), eprint 1710.06168.
  * Pardo and Spergel (2020) K. Pardo and D. N. Spergel, Phys. Rev. Lett. 125, 211101 (2020), eprint 2007.00555.
  * Ade et al. (2016) P. A. R. Ade et al. (Planck), Astron. Astrophys. 594, A27 (2016), eprint 1502.01598.
  * Gopika and Desai (2021) K. Gopika and S. Desai, arXiv e-prints arXiv:2106.07294 (2021), eprint 2106.07294.
  * Allen et al. (2011) S. W. Allen, A. E. Evrard, and A. B. Mantz, Ann. Rev. Astron. Astrophys. 49, 409 (2011), eprint 1103.4829.
  * Vikhlinin et al. (2014) A. A. Vikhlinin, A. V. Kravtsov, M. L. Markevich, R. A. Sunyaev, and E. M. Churazov, Physics Uspekhi 57, 317-341 (2014).
  * Halenka and Miller (2020) V. Halenka and C. J. Miller, Phys. Rev. D 102, 084007 (2020), eprint 1807.01689.
  * Gopika and Desai (2020) K. Gopika and S. Desai, Physics of the Dark Universe 30, 100707 (2020), eprint 2006.12320.
  * Sarazin (1986) C. L. Sarazin, Reviews of Modern Physics 58, 1 (1986).
  * Lin et al. (2012) Y.-T. Lin, S. A. Stanford, P. R. M. Eisenhardt, A. Vikhlinin, B. J. Maughan, and A. Kravtsov, Astrophys. J. Lett. 745, L3 (2012), eprint 1112.1705.
  * Rahvar and Mashhoon (2014) S. Rahvar and B. Mashhoon, Phys. Rev. D 89, 104011 (2014), eprint 1401.4819.
  * Ebeling et al. (1998) H. Ebeling, A. C. Edge, H. Bohringer, S. W. Allen, C. S. Crawford, A. C. Fabian, W. Voges, and J. P. Huchra, Mon. Not. R. Astron. Soc. 301, 881 (1998), eprint astro-ph/9812394.
  * Ebeling et al. (2000) H. Ebeling, A. C. Edge, S. W. Allen, C. S. Crawford, A. C. Fabian, and J. P. Huchra, Mon. Not. R. Astron. Soc. 318, 333 (2000), eprint astro-ph/0003191.
  * Dahle (2006) H. Dahle, Astrophys. J. 653, 954 (2006), eprint astro-ph/0608480.
  * Poole et al. (2006) G. B. Poole, M. A. Fardal, A. Babul, I. G. McCarthy, T. Quinn, and J. Wadsley, Mon. Not. R. Astron. Soc. 373, 881 (2006), eprint astro-ph/0608560.
  * Hudson et al. (2010) D. S. Hudson, R. Mittal, T. H. Reiprich, P. E. J. Nulsen, H. Andernach, and C. L. Sarazin, Astron. & Astrophys. 513, A37 (2010), eprint 0911.0409.
  * Vikhlinin et al. (2006) A. Vikhlinin, A. Kravtsov, W. Forman, C. Jones, M. Markevitch, S. S. Murray, and L. Van Speybroeck, Astrophys. J. 640, 691 (2006), eprint astro-ph/0507092.
  * Gupta and Desai (2019) S. Gupta and S. Desai, Classical and Quantum Gravity 36, 105001 (2019), eprint 1811.09378.
  * Gupta and Desai (2020) S. Gupta and S. Desai, Physics of the Dark Universe 28, 100499 (2020), eprint 1909.07408.
  * Foreman-Mackey et al. (2013) D. Foreman-Mackey, D. W. Hogg, D. Lang, and J. Goodman, Publ. Astron. Soc. Pac. 125, 306 (2013), eprint 1202.3665.
  * Hinton et al. (2020) S. Hinton, C. Adams, , and C. Badger, _Samreay/chainconsumer v0.33.0_ (2020), URL https://doi.org/10.5281/zenodo.4280904.


# Inverse Solvability and Security
with Applications to Federated Learning

Tomasz Piotrowski, Matthias Frey, Renato L.G. Cavalcante, Rafail Ismayilov
This work was in part supported by the National Centre for Research and
Development of Poland (NCBR) under grant EIG CONCERT-JAPAN/04/2021 and by the
Federal Ministry of Education and Research of Germany (BMBF) under grant
01DR21009.Tomasz Piotrowski is with Nicolaus Copernicus University in Toruń,
Poland. (e-mail: tpiotrowski@umk.pl). Matthias Frey is with Technische
Universität Berlin, Germany. Renato L. G. Cavalcante and Rafail Ismayilov are
with Fraunhofer Heinrich Hertz Institute, Berlin, Germany.

###### Abstract

We introduce concepts of inverse solvability and security for a generic linear
forward model and demonstrate how they can be applied to models used in
federated learning in wireless communications  that leverage OTA (Over-the-
Air) computation. We provide examples of such models which differ in the
resulting inverse solvability and security as defined in this paper. We also
show how the large number of users participating in a given iteration of
federated learning can be leveraged to increase both solvability and security.
Finally, we discuss possible extensions of the presented concepts to the case
of meshed networks as well as to the generic nonlinear forward model.

###### Index Terms:

inverse solvability, inverse security, federated learning

##  I Introduction

The issue of security of federated learning (FL) has been increasingly
receiving attention in recent years, and many concepts are being introduced
that aim at leveraging it [1, 2, 3, 4, 5, 6, 7]. In this paper, we use ideas
from inverse problem theory to introduce concepts of inverse solvability and
security, which may be used in tandem with other security-enhancing concepts
in FL.  We are in particular interested in methods that can be applied to the
OTA FL scheme proposed in [8]. The main purpose of this scheme is to allow for
a scaling behavior of the communication complexity that is favorable in
wireless networks with a massive number of nodes. In this work, we consider a
model that specializes to the one introduced in [8]. Namely, let us consider
the server FL model

$y=\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}f(x)+\gamma,$
(1)

where $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}f\colon\mathbb{R}^{Md}\to\mathbb{R}^{s}$
represents the ($M$-clients)$\to$server transmission model (forward model)
consisting of client-side processing (e.g., compression, sparsification) and
wireless transmission channel, $d$ is the length of the parameter vector, and
$s$ is the number of receiver endpoints (e.g., number of antennas times
frequency slots). The vector $x\in\mathbb{R}^{Md}$ represents parameter
vectors across $M$ users. The random vector $\gamma$ is spatially white with
each coefficient following $\mathcal{N}(0,\sigma^{2}_{\gamma}).$ The idea of
this model is that $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}f$
represents the pre-processing at the transmitters and the way in which the
transmitted signals are superimposed in the wireless channel, and $\gamma$
represents additive channel noise. The FL model (1) is ill-posed because $s\ll
Md$ owing to the large number of transmitted parameters $d$ and client devices
$M$ taking part in a given iteration of federated learning.

We note that FL models are usually constructed using local, deterministic and
nonlinear preprocessing operators (e.g., sparsification) applied uniformly
across users, followed by linear operators modeled as random, such as the
compression operator applied locally at each user and the global wireless
transmission channel, see, e.g., [8, 9] for up-to-date wireless FL models of
such kind. That is to say, hereafter we assume the following stochastic
representation of the forward model ${}_{M}f$ in (1):

$\forall x\in\mathbb{R}^{Md}\quad\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}f(x)\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L(H(x)),$
(2)

where

$H\colon\mathbb{R}^{Md}\to\mathbb{R}^{Md}:(x_{1},\dots,x_{M})\mapsto(h(x_{1}),\dots,h(x_{M})),$

$h\colon\mathbb{R}^{d}\to\mathbb{R}^{d}$ represents the local, deterministic
and nonlinear preprocessing operator applied uniformly across all users
$m=1,\dots,M$, and $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L$ is a
random matrix taking values in $\mathbb{R}^{s\times Md}.$ In particular, this
random matrix can be interpreted as a random linear operator composed of local
and global linear operators across $M$ users. For later use, we define the
condition number of $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L$ in a
standard way as

$\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}\sigma_{1}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)/\sigma_{s}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)\geq 1,$
(3)

where $\sigma_{1}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)$ and
$\sigma_{s}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)$ are
the largest and the smallest singular values of $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L$,
respectively. In practical settings, one has $\sigma_{s}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)>0$ with
probability 1, in which case $\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)<+\infty.$

##  II Definitions of Inverse solvability and security

The idea of inverse solvability is that the difficulty of reconstructing
preprocessed signal $H(x)$ from $y$ in model (1) with decomposition (2) may
only decrease with increasing number of users $M$ at a given iteration of an
FL scheme. Similarly, the idea of inverse security is that, for any number of
users $M$, the difficulty of reconstructing $H(x)$ should be greater for an
eavesdropper than for a legitimate server. To formalize these concepts, we
recall that, in linear inverse problems, the condition number is a common
means of identifying difficult problem instances, in the sense that the
reconstruction of $H(x)$ becomes increasingly difficult as the condition
number of $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L$
increases [10].

Based on the above considerations, we define the following notions of inverse
solvability and inverse security for the random operators in (2):

######  Definition 1 (inverse solvability).

Let $M$ be the number of active clients and let
$F\colon\mathbb{N}\to[1,\infty)$ be nonincreasing. We call the (forward) FL
model at least $F$-inverse solvable if

$\forall M\in\mathbb{N}\quad\mathbb{E}[\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)]\leq
F(M).$ (4)

######  Definition 2 (inverse security).

Assume that an eavesdropper observes the following forward model:

$y_{E}=\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}f_{E}(x)+\gamma_{E},$
(5)

yielding the same FL model representation for an eavesdropper as the server FL
model (1) with decomposition (2). Then, model (1) is inverse secure w.r.t. (5)
if

$\forall M\in\mathbb{N}\quad\mathbb{E}[\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L_{E})]>\mathbb{E}[\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)],$ (6)

where $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L_{E}$ is
the (random) linear operator observable by the eavesdropper.

##  III Examples of Inverse Solvability of FL models

We introduce first the following FL model which will be used below as a
starting point for a refined model given in Example 2, aimed in particular at
inversely secure recovery of an averaged stochastic gradients computed across
users. Similar FL models have been considered, e.g., in [8, 9].

######  Example 1.

Consider an FL model where the general aim is to recover the stochastic
gradients computed locally by $M$ users at a given iteration (time instant).
Namely,

$x=(g_{1}(\theta),\dots,g_{M}(\theta)),$ (7)

where $\theta$ is the vector of neural network parameters and $g_{m}(\theta)$
are stochastic gradients of the cost function computed using local training
data available to each of the $M$ users.

On a client side, the gradients to be transmitted are:

  1. 1.

sparsified using $S(x_{m})=g_{m}^{sp}(\theta)$, where
$S\colon\mathbb{R}^{d}\to\mathbb{R}^{d}$ is such that
$[g_{m}^{sp}(\theta)]_{i}=0$ if $|[g_{m}(\theta)]_{i}|<\delta$ and
$[g_{m}^{sp}(\theta)]_{i}=[g_{m}(\theta)]_{i}$ otherwise, for $i=1,2,\dots,d$
and where $\delta>0$ is the sparsification factor,

  2. 2.

compressed using a matrix $A\in\mathbb{R}^{s\times d}$, where $s<d$ and the
entries of $A$ are i.i.d. distributed as $\mathcal{N}(0,1).$

Denoting by $\sqrt{\alpha_{m}(t)}>0$ the power scaling coefficient of the
$m$-th (client) device, we thus consider the following FL model

$y=\frac{1}{M}A\sum_{m=1}^{M}\sqrt{\alpha_{m}}g_{m}^{sp}(\theta)+\gamma,\quad[A]_{i,j}\sim\mathcal{N}(0,1).$
(8)

This model corresponds to an additive white Gaussian noise (AWGN) multiple-
access channel (MAC). The signal transmitted by user $m$ is
$A\sqrt{\alpha_{m}}g_{m}^{sp}(\theta)$. These signals are added in the channel
together with an additive noise vector $\gamma$. Finally, the factor $1/M$ can
be seen as a normalization applied by the receiver in post-processing. Thus,
the resulting communication scheme is as proposed in [8].

######  Remark 1.

For model (8), we have in terms of decomposition (2):

$h=S,\quad H=(\underbrace{S,\dots,S}_{M}),$

and

$\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L=\frac{1}{M}(\sqrt{\alpha_{1}}A,\dots,\sqrt{\alpha_{M}}A),\quad[A]_{i,j}\sim\mathcal{N}(0,1).$

An assumption on distribution of power scaling coefficients
$\sqrt{\alpha_{m}(t)}>0$ is made neither in model (8) nor in any of the FL
models below.

######  Proposition 1.

The model in Example 1 is at least $F$-inverse solvable with
$F(M)\equiv\frac{\sqrt{d}+\sqrt{s}}{\sqrt{d}-\sqrt{s}}$ and
$\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)\equiv\operatorname{cond}(A).$

###### Proof.

We first note that the condition number does not depend on rescaling a given
matrix by a fixed scalar, so it is sufficient to investigate the condition
number of

$\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime}=(\sqrt{\alpha_{1}}A,\dots,\sqrt{\alpha_{M}}A).$

Let $\sigma_{i}(X)$ and $\lambda_{i}(Y)$ denote the $i$-th singular value and
eigenvalue of a given matrix $X$ and a symmetric matrix $Y$, respectively,
both organized in nonincreasing order. Then

$\sigma_{i}^{2}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})=\lambda_{i}((\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})^{t})=\\\
\lambda_{i}(\sum_{m=1}^{M}\alpha_{m}AA^{t})=\sigma_{i}^{2}(A)\sum_{m=1}^{M}\alpha_{m}.$

Let
$\alpha\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}\sum_{m=1}^{M}\alpha_{m}.$
Then $\sigma_{i}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})=\sqrt{\alpha}\sigma_{i}(A)$
for $i=1,2,\dots,s.$ Therefore,

$\forall M\in\mathbb{N}\quad\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)=\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})=\\\
\frac{\sigma_{1}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)}{\sigma_{s}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)}=\frac{\sigma_{1}(A)}{\sigma_{s}(A)}=\operatorname{cond}(A),$

hence also $\mathbb{E}[\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)]=\mathbb{E}[\operatorname{cond}(A)].$111We
note that, for $[A]_{i,j}\sim\mathcal{N}(0,1)$,
$\sigma_{s}(A)\geq\sqrt{d}-\sqrt{s}$ with probability 1, see [11] and, in
particular, [12, Eq. (2.3)]. Moreover, from Fact 1 in the Appendix one has
that

$\sqrt{d}-\sqrt{s}\leq\mathbb{E}\sigma_{s}(A)\leq\mathbb{E}\sigma_{1}(A)\leq\sqrt{d}+\sqrt{s},$

which yields the desired assertion. $\blacksquare$

In the signal processing literature, it is frequently assumed that the server
distributes the compression matrix $A$ in model (8) to the users, see, e.g.,
[8, 9] and references therein. However, in many FL scenarios, a large number
$M$ of users participates in a given learning iteration. We can exploit this
fact to provide a setup where users generate compression matrices locally,
which further enhances both security and solvability. To this end, we first
note that, from a practical perspective, it is natural to require that

$\mathbb{E}[y]=\frac{1}{M}\sum_{m=1}^{M}g_{m}^{sp}(\theta).$ (9)

Thus, to satisfy (9), we modify model (8) as follows:

######  Example 2.

Consider the following FL model:

$y=\frac{1}{M}\sum_{m=1}^{M}B_{m}\sqrt{\alpha_{m}}g_{m}^{sp}(\theta)+\gamma,\\\
[B_{m}]_{i,j}\sim\mathcal{N}(0,1/\alpha_{m}),\ i\neq j,\\\
[B_{m}]_{i,i}\sim\mathcal{N}(1/\sqrt{\alpha_{m}},1/\alpha_{m}),$ (10)

where the entries of $B_{m}$ are independent. Compared with model (8), this
amounts to redefining $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L$ as

$\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L=\frac{1}{M}(\sqrt{\alpha_{1}}B_{1},\dots,\sqrt{\alpha_{M}}B_{M}).$
(11)

This is a modification of the communication scheme from Example 1 in which the
transmitters do not employ a common compression matrix $A$, but rather each
transmitter $m$ uses its own compression matrix $B_{m}$.´

######  Proposition 2.

The model in Example 2 satisfies condition (9), and is at least $F$-inverse
solvable for

$F(M)=\frac{\sqrt{Md}+\sqrt{s}+\sqrt{M}}{\sqrt{Md}-\sqrt{s}-\sqrt{M}}$ (12)

if $(\sqrt{d}-1)^{2}>s.$ In such a case, $F(M)$ is monotonically decreasing.
Moreover,

$\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)=\operatorname{cond}(C_{1}+I_{s\times
d},\dots,C_{M}+I_{s\times d}),\\\ [C_{m}]_{i,j}\sim\mathcal{N}(0,1),$ (13)

where $I_{s\times d}$ is the matrix with ones on the main diagonal and zeros
elsewhere, and the entries of $C_{m}$ are independent.

###### Proof.

We first note that

$\mathbb{E}[y]=\frac{1}{M}\sum_{m=1}^{M}\mathbb{E}[B_{m}\sqrt{\alpha_{m}}g_{m}^{sp}(\theta)]$

due to linearity of the expectation operator and $\mathbb{E}[\gamma]=0.$ Let
$m\in\\{1,2,\dots,M\\}$ and $i\in\\{1,2,\dots,s\\}.$ Then

$\mathbb{E}[B_{m}\sqrt{\alpha_{m}}g_{m}^{sp}(\theta)]_{i}=\mathbb{E}[\sqrt{\alpha_{m}}\sum_{j=1}^{d}[B_{m}]_{i,j}[g_{m}^{sp}(\theta)]_{j}]=\\\
\sqrt{\alpha_{m}}\sum_{j=1}^{d}\mathbb{E}[[B_{m}]_{i,j}][g_{m}^{sp}(\theta)]_{j}=[g_{m}^{sp}(\theta)]_{i}.$

This implies (9) for model (10). Let us now consider the condition number of
$\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L$ in
(11), or, equivalently, the condition number of

$\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime}=(\sqrt{\alpha_{1}}B_{1},\dots,\sqrt{\alpha_{M}}B_{M}).$

The random matrices $B_{m}$ can be decomposed into
$B_{m}=A_{m}+\frac{1}{\sqrt{\alpha_{m}}}I_{s\times d}$, where
$[A_{m}]_{i,j}\sim\mathcal{N}(0,1/\alpha_{m}).$ Thus, the random matrices
$\sqrt{\alpha_{m}}B_{m}$ satisfy
$\sqrt{\alpha_{m}}[B_{m}]_{i,j}\sim\mathcal{N}(0,1)$ for $i\neq j$ and
$\sqrt{\alpha_{m}}[B_{m}]_{i,i}\sim\mathcal{N}(1,1)$ for $i=1,\dots,s$ and
$j=1,\dots,d.$ Therefore, we can rewrite $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime}$
as

$\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime}=(C_{1}+I_{s\times
d},\dots,C_{M}+I_{s\times d}),\quad[C_{m}]_{i,j}\sim\mathcal{N}(0,1).$ (14)

Now, from [13, Theorem 3.3.16], we conclude that

$\sigma_{1}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})\leq\sigma_{1}(C_{1},\dots,C_{M})+\sigma_{1}(I_{s\times
d},\dots,I_{s\times d})=\\\ \sigma_{1}(C_{1},\dots,C_{M})+\sqrt{M}.$ (15)

As a result, from Fact 1 in the Appendix one has that

$\mathbb{E}\sigma_{1}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})\leq\sqrt{Md}+\sqrt{s}+\sqrt{M}.$
(16)

Similarly, from [13, Theorem 3.3.16], we also have that

$\sigma_{s}(C_{1},\dots,C_{M})\leq\sigma_{s}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})+\sqrt{M}.$
(17)

Since $\mathbb{E}\sigma_{s}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})\leq\mathbb{E}\sigma_{1}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})$,
from (17) and Fact 1 in the Appendix one also has that

$\sqrt{Md}-\sqrt{s}-\sqrt{M}\leq\mathbb{E}\sigma_{s}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime}).$
(18)

It is simple to verify that $\sqrt{Md}-\sqrt{s}-\sqrt{M}>0$ for
$M\in\mathbb{N}$ if $(\sqrt{d}-1)^{2}>s.$ Merging (16) and (18) we conclude
that

$\mathbb{E}[\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)]=\mathbb{E}[\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L^{\prime})]\leq\\\
\frac{\sqrt{Md}+\sqrt{s}+\sqrt{M}}{\sqrt{Md}-\sqrt{s}-\sqrt{M}}.$ (19)

The monotonic decrease of $F(M)$ in (12) can now be verified by extending it
to the interval $M\in[1,+\infty)$ and checking that its derivative is
negative, e.g., using Mathematica command Reduce[F’[M]<0 && d>=s>=1 &&
(Sqrt[d]-1)2>s && M>=1,M]. $\blacksquare$

We close this section with the following proposition, which establishes the
probability that $y$ in (10) is an approximation of $\mathbb{E}[y]$ in (9)
with a predefined accuracy.

######  Proposition 3.

Consider the FL model in Example 2 and let

$z\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}\frac{1}{M^{2}}\sum_{m=1}^{M}\|g_{m}^{sp}(\theta)\|^{2}+\sigma^{2}_{\gamma}>0.$

Then

$\mathbb{P}\Big{(}\frac{\|y-\mathbb{E}[y]\|^{2}}{d}\geq\epsilon\Big{)}=1-F_{\chi^{2}_{d}}\Big{(}\frac{d\epsilon}{z}\Big{)}\leq
e^{-\beta d},$ (20)

where $F_{\chi^{2}_{d}}\colon\mathbb{R}\to[0,1]$ is the cumulative
distribution function of $\chi^{2}_{d}$-distributed random variable and
$\beta\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}\min\\{\mu,\mu^{2}\\}$
with
$\mu\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}\frac{1}{4}(\frac{\epsilon}{z}-1).$

###### Proof.

Let $i\in\\{1,\dots,s\\}.$ We note first that $[y]_{i}$ is Gaussian as a
linear combination of Gaussian random variables and that, using assumptions of
model (10) and similar arguments as in the proof of Proposition 2, one obtains
$\operatorname{Var}([y]_{i})=z.$ Therefore, the random vector
$(y-\mathbb{E}[y])/\sqrt{z}$ contains $d$ i.i.d. components following
$\mathcal{N}(0,1)$ and consequently,

$Y\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}\frac{\|y-\mathbb{E}[y]\|^{2}}{z}\sim\chi^{2}_{d}.$

Thus,

$\mathbb{P}\Big{(}\frac{\|y-\mathbb{E}[y]\|^{2}}{d}\geq\epsilon\Big{)}=\mathbb{P}\Big{(}\frac{\|y-\mathbb{E}[y]\|^{2}}{z}\geq\frac{d\epsilon}{z}\Big{)}=\\\
\mathbb{P}\Big{(}Y-d\geq d(\frac{\epsilon}{z}-1)\Big{)}.$

The equality in (20) is established using the first equality above. To
complete the proof, we set now
$\delta\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}d(\frac{\epsilon}{z}-1)$
in Fact 2 in the Appendix to obtain the inequality in (20). $\blacksquare$

A numerical evaluation of (20) shows that, on average, for sample values
$d=1000$ and $s=d/2$, $M>3000$ is sufficient to obtain sufficiently small
approximation error in (20). This is in line with common benchmarks for FL
considered in literature (see [14] for recent examples) which use, e.g.,
$M=3400$ or $M=342477.$

##  IV Examples of Inverse Security of FL models

Let us return to the FL model of Example 1 and assume that the eavesdropper
has access to the following model, recognized as a perturbation of model (8):

######  Example 3.

Let us modify the FL model in (8) by considering the following model available
to the eavesdropper:

$y_{E}=\frac{1}{M}\sum_{m=1}^{M}\sqrt{\alpha_{m}}H_{E_{m}}Ag_{m}^{sp}(\theta)+\gamma,\quad[A]_{i,j}\sim\mathcal{N}(0,1),$
(21)

where matrices $H_{E_{m}}\in\mathbb{R}^{s\times s}$ induce mismatch on users’
data incurred at eavesdropper’s site compared with server’s site. No
assumption on distribution of $H_{E_{m}}$ is made in model (21).  This is a
generalization of Example 1 in which channel fading is also accounted for by
introducing the fading matrices $H_{E_{m}}$.

######  Proposition 4.

The model in Example 1 is inverse secure if the model in Example 3 is
available to the eavesdropper and

$\mathbb{E}[\operatorname{cond}(H_{E_{1}},\dots,H_{E_{M}})]>\frac{\max_{m}\sqrt{\alpha_{m}}}{\min_{m}\sqrt{\alpha_{m}}}(\mathbb{E}[\operatorname{cond}(A)])^{2}.$
(22)

###### Proof.

Following the notation in Proposition 1, one has in the case of model (21):

$\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L_{E}^{\prime}=(\sqrt{\alpha_{1}}H_{E_{1}}A,\dots,\sqrt{\alpha_{M}}H_{E_{M}}A)=P_{M}Q,$

where
$P_{M}\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}(H_{E_{1}}A,\dots,H_{E_{M}}A)\in\mathbb{R}^{s\times
Md}$ and
$Q\stackrel{{\scriptstyle\mathclap{\tiny\mbox{df.}}}}{{=}}\operatorname{diag}({\sqrt{\alpha_{1}}I_{d\times
d},\dots,\sqrt{\alpha_{M}}I_{d\times d}})\in\mathbb{R}^{Md\times Md}$ is a
block-diagonal matrix. Therefore, from Fact 3 in the Appendix, one has

$\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L_{E}^{\prime})\geq\operatorname{cond}(P_{M})\operatorname{cond}(Q)^{-1}=\\\
\operatorname{cond}(P_{M})\frac{\min_{m}\sqrt{\alpha_{m}}}{\max_{m}\sqrt{\alpha_{m}}}.$
(23)

We now reapply Fact 3 in the Appendix to obtain

$\operatorname{cond}(P_{M})\geq\\\
\operatorname{cond}(H_{E_{1}},\dots,H_{E_{M}})(1/\operatorname{cond}(\operatorname{diag}({A,\dots,A})))=\\\
\operatorname{cond}(H_{E_{1}},\dots,H_{E_{M}})(1/\operatorname{cond}(A)).$
(24)

By merging inequalities (23) and (24), and using Proposition 1 we obtain (22).
$\blacksquare$

######  Remark 2.

The main problem with an inequality such as the one in (6) for the model in
Example 1 is that $\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)\equiv\operatorname{cond}(A)$
[see Proposition 1]. As an example, if
$[H_{E_{m}}]_{i,j}\sim\mathcal{N}(0,1)$, similarly as in the proof of
Proposition 2, it can be verified using Fact 1 in the Appendix that

$\mathbb{E}[\operatorname{cond}(H_{E_{1}},\dots,H_{E_{M}})]\leq\frac{\sqrt{Md}+\sqrt{s}}{\sqrt{Md}-\sqrt{s}}.$

In such a case, there exists $M_{0}\in\mathbb{N}$ such that for $M>M_{0}$ the
condition (22) does not hold, and the inverse security is not guaranteed by
Proposition 4 even if $\sqrt{\alpha_{m}}\equiv\sqrt{\alpha}$ across $M$ users.

######  Example 4.

Let us modify the FL model in (10) by considering the following model
available to the eavesdropper:

$y_{E}=\frac{1}{M}\sum_{m=1}^{M}H_{E_{m}}B_{m}\sqrt{\alpha_{m}}g_{m}^{sp}(\theta)+\gamma,\\\
[B_{m}]_{i,j}\sim\mathcal{N}(0,1/\alpha_{m}),\ i\neq j,\\\
[B_{m}]_{i,i}\sim\mathcal{N}(1/\sqrt{\alpha_{m}},1/\alpha_{m}),$ (25)

where matrices $H_{E_{m}}\in\mathbb{R}^{s\times s}$ induce mismatch on users’
data incurred at eavesdropper’s site compared with server’s site. Similarly as
above, no assumption on distribution of $H_{E_{m}}$ is made in model (25).

######  Proposition 5.

Model in Example 2 is inverse secure if model in Example 4 is available to the
eavesdropper and condition (27) provided below holds.

###### Proof.

Following notation of Proposition 2, one has in the case of model (25)

$\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L_{E}=\frac{1}{M}(\sqrt{\alpha_{1}}H_{E_{1}}B_{1},\dots,\sqrt{\alpha_{M}}H_{E_{M}}B_{M}),$

and consequently (see derivation of (14)),

$\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L_{E}^{\prime}=(H_{E_{1}}(C_{1}+I_{s\times
d}),\dots,H_{E_{M}}(C_{M}+I_{s\times d})),\\\
[C_{m}]_{i,j}\sim\mathcal{N}(0,1).$ (26)

Therefore, using Proposition 2, we conclude that a sufficient condition for
(6) is in this case

$\forall M\in\mathbb{N}\\\
\mathbb{E}[\operatorname{cond}(H_{E_{1}}(C_{1}+I_{s\times
d}),\dots,H_{E_{M}}(C_{M}+I_{s\times d}))]>\\\
\mathbb{E}[\operatorname{cond}(C_{1}+I_{s\times d},\dots,C_{M}+I_{s\times
d})],\\\ [C_{m}]_{i,j}\sim\mathcal{N}(0,1).\ \blacksquare$ (27)

We note that, unlike condition (22) in Proposition 4, condition (27) in
Proposition 5 does not depend on transmit power of users and that both sides
of inequality (27) depend on the number of users $M.$ We also note that, for
$[C_{m}]_{i,j}\sim\mathcal{N}(0,1)$, if also
$[H_{E_{m}}]_{i,j}\sim\mathcal{N}(0,1)$, then matrices $H_{E_{m}}$ induce
correlation across columns of $H_{E_{m}}(C_{m}+I_{s\times d})$ compared with
$C_{m}+I_{s\times d}.$ This observation is verified numerically in Fig.1.

Figure 1: Sample $\operatorname{cond}(C_{1}+I_{s\times
d},\dots,C_{M}+I_{s\times d})$ and
$\operatorname{cond}(H_{E_{1}}(C_{1}+I_{s\times
d}),\dots,H_{E_{M}}(C_{M}+I_{s\times d}))$ as a function of number of users
$M$ for standard i.i.d. Gaussian $C_{m}$ and $H_{E_{m}}$, with $d=1000$ and
$s=250$, averaged over 10 independent trials.

##  V Conclusion

We provided two related concepts of inverse solvability and security of FL
models which may help in design of FL models. The provided concepts are easily
extendable, for example, one may strengthen them by removing expectation
operators in inequalities (4) and (6). Similarly, an inverse security for a
scenario different from the one considered in Definition 2 may be introduced.
We provide an example below for meshed networks.

######  Remark 3 (inverse security of meshed networks).

Assume that an eavesdropper is able to access data sent by at most $M^{\star}$
users. Then, a sufficient condition of inverse security for such networks is
to require

$\forall
M\in\\{1,\dots,M^{\star}\\}\quad\mathbb{E}[\operatorname{cond}(\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{L}}_{\mathmakebox[0pt][r]{M}}}L)]>\mathbb{E}[\operatorname{cond}(_{(M^{\star}+1)}L)].$
(28)

Finally, we should also note that, as a future work, one may envision
providing definitions similar to those presented in this paper for a more
general forward model $\mathchoice{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
6.04332pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}{\kern
4.31665pt\mathop{}\mathopen{\vphantom{f}}_{\mathmakebox[0pt][r]{M}}}f$, e.g.,
with nonlinear operator $H$ varying across users. Such approach would allow to
avoid the assumptions on the forward model such as decomposition (2) used in
this paper.

######  Fact 1 ([11], see also Theorem 2.6 in [12]).

Let $X\in\mathbb{R}^{N\times n}$ be a random matrix with entries of $X$ i.i.d.
distributed as $\mathcal{N}(0,1).$ Then

$\sqrt{N}-\sqrt{n}\leq\mathbb{E}[\sigma_{\min}(X)]\leq\mathbb{E}[\sigma_{\max}(X)]\leq\sqrt{N}+\sqrt{n},$
(29)

where $\sigma_{\min}(X)$ and $\sigma_{\max}(X)$ are the smallest and the
largest singular values of $X$, respectively.

######  Fact 2 (Corollary from Lemma 1, eq. (4.1) in [15]).

Let $Y\sim\chi^{2}_{d}$ and let $\delta>0.$ Then

$\mathbb{P}(Y-d\geq\delta)\leq\exp\Big{(}-\min\\{\frac{\delta}{4},\frac{\delta^{2}}{16d}\\}\Big{)}.$
(30)

######  Fact 3 (Corollary from Theorem 3.3.16 in [13]).

Let $P\in\mathbb{R}^{m\times n}$ and $Q\in\mathbb{R}^{n\times s}$ be of full
rank and such that $m\leq n\leq s.$ Then

$\operatorname{cond}(PQ)\geq\operatorname{cond}(P)\operatorname{cond}(Q)^{-1}.$
(31)

## References

  * [1] Solmaz Niknam, Harpreet S. Dhillon, and Jeffrey H. Reed,  “Federated learning for wireless communications: Motivation, opportunities, and challenges,”  IEEE Communications Magazine, vol. 58, no. 6, pp. 46–51, 2020.
  * [2] Kang Wei, Jun Li, Ming Ding, Chuan Ma, Howard H. Yang, Farhad Farokhi, Shi Jin, Tony Q. S. Quek, and H. Vincent Poor,  “Federated learning with differential privacy: Algorithms and performance analysis,”  IEEE Transactions on Information Forensics and Security, vol. 15, pp. 3454–3469, 2020.
  * [3] Chuan Ma, Jun Li, Ming Ding, Howard H. Yang, Feng Shu, Tony Q. S. Quek, and H. Vincent Poor,  “On safeguarding privacy and security in the framework of federated learning,”  IEEE Network, vol. 34, no. 4, pp. 242–248, 2020.
  * [4] Matthias Frey, Igor Bjelaković, and Sławomir Stańczak,  “Towards secure over-the-air computation,”  Submitted to IEEE Transactions on Information Theory, 2022,  Preprint available at arXiv:2001.03174.
  * [5] Matthias Frey, Igor Bjelaković, and Sławomir Stańczak,  “Towards secure over-the-air computation,”  in 2021 IEEE International Symposium on Information Theory (ISIT). IEEE, 2021, pp. 700–705.
  * [6] Yi Liu, Jialiang Peng, Jiawen Kang, Abdullah M. Iliyasu, Dusit Niyato, and Ahmed A. Abd El-Latif,  “A secure federated learning framework for 5G networks,”  IEEE Wireless Communications, vol. 27, no. 4, pp. 24–31, 2020.
  * [7] Changjie Hu, Quanzhong Li, Qi Zhang, and Jiayin Qin,  “Secure transceiver design and power control for over-the-air computation networks,”  IEEE Communications Letters, vol. 26, no. 7, pp. 1509–1513, 2022\.
  * [8] Mohammad Mohammadi Amiri and Deniz Gündüz,  “Machine learning at the wireless edge: Distributed stochastic gradient descent over-the-air,”  IEEE Transactions on Signal Processing, vol. 68, pp. 2155–2169, 2020.
  * [9] Mohammad Mohammadi Amiri and Deniz Gündüz,  “Federated learning over wireless fading channels,”  IEEE Transactions on Wireless Communications, vol. 19, no. 5, pp. 3546–3557, 2020.
  * [10] G. H. Golub and C. F. Van Loan,  Matrix Computations,  The Johns Hopkins University Press, Baltimore, 1996.
  * [11] K. R. Davidson and S. J. Szarek,  “Chapter 8 local operator theory, random matrices and Banach spaces,”  vol. 1 of Handbook of the Geometry of Banach Spaces, pp. 317 – 366\. Elsevier Science B.V., 2001.
  * [12] M. Rudelson and R. Vershynin,  “Non-asymptotic theory of random matrices: extreme singular values,”  in Proceedings of the International Congress of Mathematicians, 2010\.
  * [13] R. A. Horn and C. R. Johnson,  Topics in Matrix Analysis,  Cambridge University Press, New York, 1991.
  * [14] Wei-Ning Chen, Christopher A. Choquette-Choo, Peter Kairouz, and Ananda Theertha Suresh,  “The fundamental price of secure aggregation in differentially private federated learning,” 2022,  arXiv.
  * [15] B. Laurent and P. Massart,  “Adaptive estimation of a quadratic functional by model selection,”  The Annals of Statistics, vol. 28, no. 5, pp. 1302 – 1338, 2000\.


Accepted manuscript IEEE/EMBS Neural Engineering (NER) 2021

# Inter-subject Deep Transfer Learning for Motor Imagery EEG Decoding

Xiaoxi Wei1, Pablo Ortega1 and A. Aldo Faisal1,2 Brain & Behaviour Lab: 1
Department of Computing,2 Department of Bioengineering, Imperial College
London. Correspondence<EMAIL_ADDRESS>
###### Abstract

Convolutional neural networks (CNNs) have become a powerful technique to
decode EEG and have become the benchmark for motor imagery EEG Brain-Computer-
Interface (BCI) decoding. However, it is still challenging to train CNNs on
multiple subjects’ EEG without decreasing individual performance. This is
known as the negative transfer problem, i.e. learning from dissimilar
distributions causes CNNs to misrepresent each of them instead of learning a
richer representation. As a result, CNNs cannot directly use multiple
subjects’ EEG to enhance model performance directly. To address this problem,
we extend deep transfer learning techniques to the EEG multi-subject training
case. We propose a multi-branch deep transfer network, the Separate-Common-
Separate Network (SCSN) based on splitting the network’s feature extractors
for individual subjects. We also explore the possibility of applying Maximum-
mean discrepancy (MMD) to the SCSN (SCSN-MMD) to better align distributions of
features from individual feature extractors. The proposed network is evaluated
on the BCI Competition IV 2a dataset (BCICIV2a dataset) and our online
recorded dataset. Results show that the proposed SCSN (81.8%, 53.2%) and SCSN-
MMD (81.8%, 54.8%) outperformed the benchmark CNN (73.4%, 48.8%) on both
datasets using multiple subjects. Our proposed networks show the potential to
utilise larger multi-subject datasets to train an EEG decoder without being
influenced by negative transfer.

###### Index Terms:

brain-computer-interface, EEG, multi-subject, deep learning, transfer
learning, online decoding

##  I Introduction

Since its first application in BCI in the mid-2010s [1], deep learning has
emerged as a powerful technique to learn EEG decoders without biased feature
design caused by the limited understanding of human EEG. Convolutional neural
networks (CNNs) have outperformed conventional machine learning methods and
become the benchmark for motor imagery EEG decoding [2]. However, CNNs are,
unlike conventional BCI machine learning [3], not data-efficient, so combining
multiple datasets is essential. EEG features vary across subjects, for
example, due to different sensor locations, sensor impedance, and individual
brain anatomical and functional organisation. A challenge associated with CNNs
is the negative transfer [4], i.e. learning from dissimilar EEG distributions
from different subjects causes CNNs to misrepresent each of them instead of
learning a richer representation.

As a result, unlike other machine learning fields, CNNs cannot directly
utilise multiple subjects’ EEG to enhance model performance and make up for
the data shortage. Additionally, BCI systems usually require a long training
period for new users to get a personalised and stabilised decoder, which
disrupts the use of BCI in real-world. To address this problem, we combine
deep transfer learning with EEG-BCI decoding.

Transfer learning has emerged as a solution to combine data from domains with
different distributions [5]. Transfer learning encompasses algorithms which
aim to transfer the representations and knowledge from source domains to a
target domain. There have been some promising strategies for deep transfer
learning in fields like computer vision and natural language processing. For
example, fine-tuning is often used to transfer model representations from one
task or dataset to another [6], i.e. networks are pre-trained on one domain
and then retrained on another domain. Another strategy is to train a shared
network using multiple datasets but split deep network layers for different
datasets [7]. In addition, deep domain adaptation [8] is a technique which
directly projects source and target distributions into a common distribution
with CNNs, and thus data from source domains can benefit the model performance
on the target domain.

Figure 1: Illustration of Deep Learning and Deep Transfer Learning in BCI.

Few studies have attempted to adapt the ideas above in EEG decoding. In [9], a
deep transfer network was proposed to transfer the knowledge representation
from images to EEG. Following this study, an adversarial network was added
into the network to perform better domain adaptation [10]. Another study used
a naive idea of pre-training and fine-tuning [11] for motor imagery decoding
models. A deep domain adaptation network was used to transfer the data
representation from one subject to another [12], in which deeper layers of the
network are separated for the two individuals. However, this domain adaptation
network only transfers one subject’s distribution to another subject (one-to-
one). This limits the amount of data that can be used to train a deep learning
model. It is still challenging in the literature for the above deep transfer
learning methods to transfer across multiple subjects.

In this study, we designed a multi-subject deep transfer learning network and
tested it on the BCI Competition IV 2a Dataset (BCICIV2a dataset) [13] and our
own live, online dataset.

##  II Methods

To begin with, we use the benchmark motor imagery decoder [2] as a baseline.
The CNN decoder consists of a temporal layer, a spatial layer, a mean pooling
layer and a classification layer. The temporal layer extracts time relevant
information. The spatial layer performs a spatial convolution across EEG
channels. The mean pooling layer averages across features, and finally, the
classification layer uses the averaged features to predict a classification
label.

We design a Separate-Common-Separate Network (SCSN) by separating feature
extractor of the baseline CNN for individuals, so each of the subjects has
their own temporal layer, spatial layer and mean pooling layer, as shown in
Figure 2. Each network branch extracts subject-specific features, which could
avoid negative transfer in the common feature extractor of the baseline CNN.
After that, there are three fully connected layers to extract common features
of all subjects in the feature space. Unlike the baseline CNN, SCSN extracts
common features in the feature space, which could avoid learning from the
noise of raw data. Differences in individual’s EEG are explained by many
aspects, e.g. sensor positions during recording, sensor impedance, and
individual brain functionalities. Intuitively, splitting shallow feature
extractors could deal with features more related to raw data, i.e. sensor
positions and impedance. EEG features which reveal brain functionalities could
present in deep layers of the network after several layers of feature
extraction. In light of this, the SCSN separate the feature extractors again
in deeper layers before the classification layer to handle differences in
individual brain functionality.

Figure 2: Architecture of the Separate-Common-Separate Network (SCSN).

$\mbox{MMD}_{i}^{2}(X_{si},X_{t})=\\\
||\frac{1}{{X_{s}}_{i}}\sum_{{x_{s}}_{i}\in{X_{s}}_{i}}\Phi({x_{s}}_{i})-\frac{1}{X_{t}}\sum_{x_{t}\in
X_{t}}\Phi(x_{t})||.$

Maximum-mean discrepancy (MMD) [14] is a metric which measures the distance
between two datasets in kernel space. For the Separate-Common-Separate Network
with MMD (SCSN-MMD), we compute MMD between the target and each source subject
in the separate deep feature extractors and add it into the loss function.
Equation II shows the computation of MMD. According to the literature [15], we
use an RBF kernel with the mean L2 distance of data as the variance.
$L_{c}+\lambda\sum_{i=1}^{N}\mbox{MMD}_{i}^{2}(X_{si},X_{t})$ shows the loss
function, where $L_{c}$ represents the classification loss, N represents the
number of source subjects and $\lambda$ represents the balancing factor
between MMD loss and the classification loss. In addition, SCSN has three
fully connected layers for individual subjects before the classification
layer. We compute MMD for each of the three layers across subjects. The
$\mbox{MMD}_{i}^{2}$ loss is a weighted average MMD of each of the three
layers. Their averaging weights are $\frac{1}{6},\frac{1}{3},\frac{1}{2}$
respectively; this is to increase the significance of deeper layers. We also
match samples with the same label when we compute MMD. This ensures the MMD
loss represents the distributional distance between data with the same label
from different subjects. By projecting the target subject (the new BCI user)
and source subjects (supplementary subjects) into a similar domain, the MMD
constraint could benefit the model performance of the target subject with more
data.

##  III Evaluation

###  III-A Dataset

The BCI Competition IV 2a Dataset (BCICIV2a dataset [13]) dataset is
frequently used in the literature for motor imagery decoding. The dataset
consists of nine subjects. For each subject, two sessions were recorded, and
each session consists of 288 four-second trials of four motor imageries (left
hand, right hand, both feet and tongue). There are 22 channels covering the
scalp. Five of the nine subjects with highest data quality (subject A01, A03,
A07, A08, A09) are selected from the dataset according to published results of
the competition [13].

Additionally, we have collected our own online recorded dataset using the
Cybathlon 2020 BCI game (https://cybathlon.ethz.ch/en/event/disciplines/bci)
with the same hardware set up as our previous work in Cybathlon 2016 [16]. The
experimental procedures involving human subjects described in this paper were
approved by the Science Engineering Technology Research Ethics Committee of
Imperial College London. In the Cybathlon BCI game, subjects attempt to
control a virtual wheelchair on a winding race track to go straight, turn
left, turn right and turn on the headlight. Our dataset consists of five
offline sessions from five right-handed subjects (male), one offline session
and one online session from our target subject (right-handed, male), i.e. the
pilot. Each session contains 300 5-second trials of four motor imageries
(relaxing, left hand, right hand and both feet) of 64 channels covering the
whole scalp. In offline sessions, the virtual wheelchair in the Cybathlon BCI
game automatically ran on the winding track and subjects were asked to imagine
the movements according to actions of the wheelchair. In the online session,
the pilot used the baseline CNN model trained on his offline session to
perform real-time control to the virtual wheelchair.

Figure 3 shows topographic EEG maps of the target subject in our dataset. The
topographies are computed with the averaged power of each channel across EEG
trials. The four columns stand for alpha&beta band (8Hz - 30Hz) topographic
maps (in dB) of the rest, left-hand (LH), right-hand (RH) and feet motor
imageries of the offline session and the online pilot session.

Figure 3: Topographies of the pilot subject in our dataset.

###  III-B Data processing pipeline

Both the BCICIV2a dataset and online dataset use the same preprocessing
pipeline. First, we perform a 50Hz notch filter and a bandpass filter between
1Hz and 100Hz. We then crop each trial into 2-second trials with an overlap of
1.9 seconds to better fit the real-time setup. Although using shorter trials
may decrease decoding accuracies, it benefits the transition time when
subjects attempt to switch from one command to another and the decoding time
of the algorithm.

After pre-processing, datasets are split into training sets, validation sets
and test sets. For the BCICIV2a dataset, we simulate the calibration period in
real-world BCI by including the target subject’s first 120 trials of the
second session into the training set. The training set consists of the first
session of all five subjects as well as the target subject’s first 120 trials
of the second session, so we train a respective model for each subject in our
experiments. The validation set contains the 120th to the 144th trials of the
target subject’s second session. The last 144 trials in the second session
form the test set. Similarly, for our online dataset, the training set of the
transfer learning models consists of the offline session of all six subjects
as well as the first 100 trials of the online session of the pilot subject.
The validation set contains the 100th to the 160th trials. The last 140 trials
in the second session form the test set. We select 11 channels covering the
motor cortex in our online dataset to reduce online decoding delay. These
channels are FC1, FC2, C3, C4, CP5, CP1, CP2, CP6, P3, Pz and P4.

##  IV Results

We tested the baseline CNN, SCSN and SCSN-MMD on both BCICIV2a dataset and our
online recorded dataset. There are three main observations from the
experiment.

First, the baseline CNN encountered a significant decrease in accuracy when it
is trained on multiple subjects. The blue box in Figure 4 is the accuracy of
baseline CNN trained on the target subject’s EEG (session one and the
calibration trials of session two). The orange box in Figure 4 is the accuracy
of baseline CNN trained on the whole training set described in the Evaluation
section (Note that each subject has a respective model). On the BCICIV2a
dataset, we can observe that the average classification accuracy of the
baseline CNN dropped significantly from 82% to 73.4% with data from multiple
subjects. Similarly, for the online pilot, the accuracy decreased from 54.7%
to 48.8%. The results support our assumption that the baseline CNN encounters
negative transfer.

The main finding of this study is that the multi-branch SCSN and SCSN-MMD
outperformed the baseline CNN on the training set with multiple subjects,
shown in figure 4. Both SCSN and SCSN-MMD used a batch size of 30 for each
branch, so $30N$ trials are forwarded to the network in one training step,
where $N$ is the number of subjects. The target subject has more samples than
other source subjects because of the additional calibration period, so we
randomly duplicate samples, with a balanced size for labels, in each training
set of the source subjects. This is to balance the learning procedure of
different branches. The average accuracies of the SCSN and SCSN-MMD reached
81.8%, which are significantly higher than the accuracy of the baseline CNN
(73.4%) using multiple subjects’ EEG. We observed similar results on our
online recorded dataset. The SCSN (53.2%) and SCSN-MMD (54.8%) outperformed
the multi-subject baseline (48.8%). Note that in the testing procedure of SCSN
and SCSN-MMD, EEG trials only go through the target branch, so the main
difference between the proposed network and the baseline is the networking
training.

Finally, comparing SCSN and SCSN-MMD, we did not observe an obvious increase
in decoding accuracies by adding MMD constraints to our multi-subject
networks. The performance of the SCSN-MMD maintained 81.8% on the BCICIV2a
dataset and increased, slightly, by 1.6% on the online dataset compared with
the SCSN. All experiments are done with an RTX 2080Ti GPU. Approximate
training time for the baseline and SCSN on multiple subjects are both 700 to
800 seconds. The training time of SCSN-MMD is around 1000 seconds. Single
trial testing time for all networks are within 0.1 second.

Figure 4: Accuracies of models on BCICIV2a Dataset. The accuracies in black
represent the accuracy of each subject, and the accuracies in white show the
mean accuracy of different subjects.

##  V Discussion

The SCSN and SCSN-MMD outperformed the baseline on multiple subjects’ data
(Figure 4). A significant observation from the results is that, for inter-
subject deep transfer learning, it is essential to have an individual feature
extractor for each subject, and this is the main differences between the
proposed networks and the baseline CNN. There are many differences in features
of EEG, for example, different channel positions and impedance. We suggest
that these inter-subject differences are closely related to raw EEG data and
should to be processed in shallow layers of the decoder. Besides, brain
functionality of different subjects could be different performing the same
motor imagery [17]. This could be revealed in deeper layers of a decoder after
feature extraction in shallow layers. Our proposed networks have separate
branches for individual subjects in both shallow and deep layers. It considers
these differences in individual feature extractors and thus outperform the
multi-subject baseline.

Although the SCSN-MMD also tackle the negative transfer problem, we did not
observe significant benefits of applying constraints in our multi-branch
networks. The training time of SCSN-MMD is around 20% longer than SCSN due to
additional calculation of MMD. Performance of the SCSN and SCSN-MMD is not
significantly different. Compared with [12], in which MMD is used to transfer
from one subject to the other, our study extends the scenario to multiple
subjects and compute MMD between the target and each source set, which
introduces extra complexity. This may be one reason why applying MMD did not
benefit the transfer. To improve the MMD calculation, we could apply multi-
kernel methods [8] in future studies, where $\Phi$ in the Equation II could be
extended to multiple Gaussian kernels instead of one to achieve a better
projection of distributions.

One limitation of this study is that, although the proposed networks
outperformed the multi-subject baseline CNN, they reach similar accuracy with
single-subject baseline. In this study, we used five subjects as the source
set, which is relatively a small amount of data for deep transfer learning.
This could be one reason why the models did not outperform the single-subject
baseline. Larger datasets should be used in future studies to achieve a higher
classification accuracy - multi-branch deep transfer learning is ideal for
this purpose. Since individual feature extractors can handle different
distributions, multi-branch networks like SCSN could be potential solutions to
make use of datasets of different subjects with varying sensor sizes, sensor
locations and devices.

##  VI Conclusion

We have designed a multi-subject deep transfer learning network, the SCSN, and
tested it on both BCICIV2a dataset and our online recorded dataset. Results
show that our proposed network can tackle the negative transfer problem in the
benchmark CNN model, and have outperformed the multi-subject baseline. We
believe that our work is an essential conceptual proof towards successfully
exploiting multi-subject information in EEG decoding using multi-branch deep
transfer learning.

## References

  * [1] I. Walker, M. Deisenroth, and A. Faisal, “Deep convolutional neural networks for brain computer interface using motor imagery,” _Imperial College, Tech Report_ , p. 68, 2015.
  * [2] R. T. Schirrmeister, J. T. Springenberg, L. D. J. Fiederer, M. Glasstetter, K. Eggensperger, M. Tangermann, F. Hutter, W. Burgard, and T. Ball, “Deep learning with convolutional neural networks for eeg decoding and visualization,” _Human brain mapping_ , vol. 38, no. 11, pp. 5391–5420, 2017\.
  * [3] A. Ferrante, C. Gavriel, and A. Faisal, “Data-efficient hand motor imagery decoding in eeg-bci by using morlet wavelets & common spatial pattern algorithms,” in _IEEE/EMBS Neural Engineering (NER)_ , vol. 7. IEEE, 2015, pp. 948–951.
  * [4] S. J. Pan and Q. Yang, “A survey on transfer learning,” _IEEE Trans. on Knowledge and Data Eng._ , vol. 22, no. 10, pp. 1345–1359, 2009.
  * [5] V. Jayaram, M. Alamgir, Y. Altun, B. Scholkopf, and M. Grosse-Wentrup, “Transfer learning in brain-computer interfaces,” _IEEE Comp. Intel. Magazine_ , vol. 11, no. 1, pp. 20–31, 2016.
  * [6] J. Yosinski, J. Clune, Y. Bengio, and H. Lipson, “How transferable are features in deep neural networks?” in _Adv. in neural information processing systems (NIPS)_ , 2014, pp. 3320–3328.
  * [7] J.-T. Huang, J. Li, D. Yu, L. Deng, and Y. Gong, “Cross-language knowledge transfer using multilingual deep neural network with shared hidden layers,” in _2013 IEEE Intl. Conf. on Acoustics, Speech and Signal Processing_. IEEE, 2013, pp. 7304–7308.
  * [8] M. Long, Y. Cao, J. Wang, and M. I. Jordan, “Learning transferable features with deep adaptation networks,” _arXiv preprint arXiv:1502.02791_ , 2015\.
  * [9] C. Tan, F. Sun, and W. Zhang, “Deep transfer learning for eeg-based brain computer interface,” in _2018 IEEE Intl. Conf. on Acoustics, Speech and Signal Processing (ICASSP)_. IEEE, 2018, pp. 916–920.
  * [10] C. Tan, F. Sun, B. Fang, T. Kong, and W. Zhang, “Autoencoder-based transfer learning in brain–computer interface for rehabilitation robot,” _Intl. J. Adv. Robotic Sys._ , vol. 16, no. 2, p. 1729881419840860, 2019.
  * [11] S. Sakhavi and C. Guan, “Convolutional neural network-based transfer learning and knowledge distillation using multi-subject data in motor imagery bci,” in _8th Intl. IEEE/EMBS Conf. on Neural Engineering (NER)_. IEEE, 2017, pp. 588–591.
  * [12] W. Hang, W. Feng, R. Du, S. Liang, Y. Chen, Q. Wang, and X. Liu, “Cross-subject eeg signal recognition using deep domain adaptation network,” _IEEE Access_ , vol. 7, pp. 128 273–128 282, 2019.
  * [13] M. Tangermann, K.-R. M”̈uller, A. Aertsen, N. Birbaumer, C. Braun, C. Brunner, R. Leeb, C. Mehring, K. J. Miller, G. Mueller-Putz _et al._ , “Review of the bci competition iv,” _Frontiers in Neuroscience_ , vol. 6, p. 55, 2012.
  * [14] A. Gretton, K. Borgwardt, M. Rasch, B. Schölkopf, and A. Smola, “A kernel method for the two-sample-problem,” _Advances in neural information processing systems_ , vol. 19, pp. 513–520, 2006.
  * [15] A. Gretton, K. M. Borgwardt, M. J. Rasch, B. Sch”̈olkopf, and A. Smola, “A kernel two-sample test,” _J Machine Learning Res._ , vol. 13, no. 1, pp. 723–773, 2012.
  * [16] P. Ortega, C. Colas, and A. A. Faisal, “Compact convolutional neural networks for multi-class, personalised, closed-loop eeg-bci,” in _IEEE Biomedical Robotics and Biomechatronics (Biorob)_. IEEE, 2018, pp. 136–141.
  * [17] A. Mierau, W. Klimesch, and J. Lefebvre, “State-dependent alpha peak frequency shifts: Experimental evidence, potential mechanisms and functional implications,” _Neuroscience_ , vol. 360, pp. 146–154, 2017.