# A SEARCH FOR CO-EVOLVING ION AND NEUTRAL GAS SPECIES IN PRESTELLAR MOLECULAR
CLOUD CORES

Konstantinos Tassis11affiliation: Max-Planck Institut für Radioastronomie,
53121 Bonn, Germany , Talayeh Hezareh11affiliation: Max-Planck Institut für
Radioastronomie, 53121 Bonn, Germany , & Karen Willacy22affiliation: Jet
Propulsion Laboratory, California Institute of Technology, Pasadena, CA 91109,
USA

###### Abstract

Comparison of linewidths of spectral line profiles of ions and neutral
molecules have been recently used to estimate the strength of the magnetic
field in turbulent star-forming regions. However, the ion (HCO+) and neutral
(HCN) species used in such studies may not be necessarily co-evolving at every
scale and density and may thus not trace the same regions. Here, we use
coupled chemical/dynamical models of evolving prestellar molecular cloud cores
including non-equilibrium chemistry, with and without magnetic fields, to
study the spatial distribution of HCO+ and HCN, which have been used in
observations of spectral linewidth differences to date. In addition, we seek
new ion-neutral pairs that are good candidates for such observations because
they have similar evolution and are approximately co-spatial in our models. We
identify three such good candidate pairs: HCO+/NO, HCO+/CO, and NO+/NO.

###### Subject headings:

ISM: molecules – ISM: clouds – ISM: magnetic fields – magnetohydrodynamics
(MHD) – stars: formation – ISM: abundances

##  1\. Introduction

Observations of large velocity dispersions in the line profiles of molecules
(line widths wider than expected from thermal motions) in interstellar
molecular clouds is a signature of localized supersonic turbulence, while
observations of the Zeeman effect (e.g., Crutcher et al., 1999) and the
polarization of dust thermal emission (e.g., Hildebrand et al., 2009) indicate
that these clouds are permeated by magnetic fields. Thus interstellar
molecular clouds need to be studied in the framework of magnetohydrodynamics.

For typical molecular cloud conditions (e.g., $B\sim 10\mu G$, density $\sim
10^{3}{\rm cm^{-3}}$, and ionized fraction $\lesssim 10^{-6}$) magnetic
forces, which are experienced directly by ionized particles, are transmitted
to the bulk of the neutral gas molecules via frequent collisions. However the
coupling of the magnetic field with the gas is not perfect. Motions of the
bulk neutral gas with frequencies higher than the typical ion-neutral
collision frequency can not be transfered to the ions. This relative motion
between charged and neutral particles, known as ambipolar drift (Mestel &
Spitzer, 1956), results in the dissipation of hydromagnetic waves in small
scales (Braginskii, 1965; Kulsrud & Pearce, 1969; Zweibel & Josafatsson, 1983)
and can cause differences between the velocity dispersion spectra of the ions
and the neutrals, as was recognized by Houde et al. (2000a, b) and explored
further by Li & Houde (2008). Indeed, Mouschovias et al. (2011), in a
comprehensive study of magnetohydrodynamic waves in weakly-ionized gases,
showed that most hydromagnetic waves with wavelengths below a critical value
cannot propagate in the neutrals because they are damped rapidly by ambipolar
diffusion, but certain long-lived modes in the neutrals exist at all
wavelengths; in contrast, waves in the ions damp very rapidly at short
wavelengths.

Observationally, the narrower line profiles of ions compared to that of co-
existing neutral species was first noticed and investigated by Houde et al.
(2000a, b), who argued that the reason of this ion line narrowing effect was
due to the presence of magnetic fields in turbulent molecular clouds that are
on average not aligned with the local flows. To test this analysis, Lai et al.
(2003) observed DR21(OH) at a high spatial resolution in H13CO+ and H13CN and
confirmed that the ion-neutral pair is co-existent within their observing
scale ($6^{\prime\prime}$) and that the ion line widths are indeed narrower
than the neutral line widths. Li & Houde (2008) re-visited this magnetically
induced ion line narrowing effect. They mapped the M17 star forming region in
${\rm HCO^{+}}$ $(J=4\rightarrow 3)$ and HCN $(J=4\rightarrow 3)$ and fitted
the lower envelopes of the velocity dispersion spectra of the observed line
profiles of the ion and neutral pair to a Kolmogorov-type power-law. They
determined the ambipolar diffusion length scale using the fitting parameters,
and used it to calculate the strength of the plane-of-sky component of the
magnetic field. Hezareh et al. (2010) further tested the Li & Houde (2008)
technique by mapping DR21(OH) in the optically thin ${\rm H^{13}CO^{+}}$ and
${\rm H^{13}CN}$ $(J=4\rightarrow 3)$ lines and obtained the turbulent
ambipolar diffusion length and magnetic field strength in that source. Hezareh
et al. (2012) used interferometry maps of the ground transition of the latter
pair of molecules in a few massive dense cores in the Cygnus X region to
verify that the method is reproducible for both single-dish and interferometry
data.

One concern in interpreting these observations is that chemical
differentiation might also be responsible for differences in the spectral
appearance of different molecular species: species with different spatial
distributions would also exhibit different velocity profiles. Although HCN and
HCO+ have been shown to be approximately co-existent in the observed molecular
clouds on the scales probed (Houde et al., 2000a, b; Lai et al., 2003; Houde
et al., 2004), their chemical co-evolution has yet to be examined from a
theoretical perspective. Here, we use both magnetic and non-magnetic dynamical
models of core collapse coupled with, and self-consistently following, non-
equilibrium chemistry both in the gas phase and on grains (Tassis et al.,
2012, hereafter paper I) to examine whether the typically observed neutrals
and ions are indeed theoretically expected to be similarly distributed in a
star-forming core. We also repeat this exercise for other frequently observed
neutral and ionized molecules, and we propose neutral-ion pairs that,
according to our models, best maintain a co-spatial distribution.

##  2\. Models

Figure 1.—  Line types and colors used to denote each of the models studied.
Solid normal-thickness red line: “reference” magnetic model; solid normal-
thickness blue line: “reference” non-magnetic model. Dotted lines: “fast”
models; dashed lines: “slow” models. Brown/purple lines: magnetic/nonmagnetic
models with temperatures differing from the “reference” models. Orange/cyan
shaded areas: variation in C/O ratio. Thin/thick solid red/blue lines:
lower/higher cosmic-ray ionization rate magnetic/non-magnetic models (see text
for details).

We are interested in the co-spatiality of various observable neutral-ion pairs
as a function of density in a star-forming region. Since the behavior of the
abundance of many molecules with density is sensitive to various model
parameters, we use all Paper I models to assess the behavior of the abundance
ratios of neutral-ion pairs with density. The best candidate molecule pairs
for studies of systematic differences of the linewidths between neutrals and
ions will be those with abundance ratios that vary weakly with density for any
set of model parameters.

The physics of the chemodynamical models we use and its effect on molecular
abundances are exhaustively discussed in Paper I. Here, we briefly review the
parameters we have varied for each class of models: temperature, C/O ratio,
cosmic ray ionization rate, and a parameter controlling the time available for
chemical evolution. The latter is the mass-to-magnetic-flux ratio for magnetic
models, and the collapse delay time for non-magnetic models (an initial time
period during which chemistry evolves but the core does not evolve
dynamically, representing an early stage of support due to turbulence which
later decays).

The initial parameter values for our “reference” magnetic model are a mass-to-
magnetic-flux ratio equal to the critical value for collapse (Mouschovias &
Spitzer, 1976), a temperature of 10 K, a C/O ratio of 0.4, and the canonical
value for the cosmic ray ionization rate $\zeta=1.3\times
10^{-17}{\rm\,s^{-1}}$. The relevant parameter values of our “reference” non-
magnetic model are a collapse delay time (hereafter “delay”) of 1 Myr, with
the other parameters being identical to the “reference” magnetic model.

For magnetic models we vary the initial mass-to-magnetic-flux ratio, examining
two additional values: 1.3 times the critical value (a faster-evolving,
magnetically supercritical model), and 0.7 of the critical value (a slower,
magnetically subcritical model). For non-magnetic models we examine two
additional values of delay: zero, and 10 Myr.

For each of these six dynamical models, the carbon-to-oxygen ratio is varied
from its “reference” value by keeping the abundance of C constant and changing
that of O. Additional values of C/O ratio we examined are 1 and 1.2. This way,
we study a “matrix” of 18 combinations (9 magnetic and 9 non-magnetic) of
available times for chemical evolution and C/O ratios.

To test the effect of temperature ($T$), we have varied each of the six basic
dynamical models by changing $T$ by a factor of $\sim 1.5$ from its reference
value of 10 K and examined models with $T=7$ K and $T=15$ K. This results in
12 additional models (6 magnetic and 6 non-magnetic).

To test the cosmic ray ionization variation, we studied four additional models
(two magnetic and two non-magnetic), which have the “reference” value for the
temperature, C/O ratio, and mass-to-flux ratio or delay (for magnetic and non-
magnetic models respectively), but for which $\zeta$ is varied by a factor of
four above ($\zeta=5.2\times 10^{-17}$ $s^{-1}$) and below ($\zeta=3.3\times
10^{-18}$ $s^{-1}$) its “reference” value, covering the range of observational
estimates (e.g., McCall et al. (2003); Hezareh et al. (2008)).

Figure 1 shows a visual representation of these suite of models and is a quick
reference guide for the line colors and types we have used to depict each
model.

##  3\. Results

###  3.1. HCO+/HCN

Figure 2.—  Left column: Evolution of central abundance ratio of HCO+ and HCN
as a function of central number density. Middle and right columns: radial
profiles of the HCO+ / HCN abundance ratio, for two different evolutionary
stages as quantified by the central number density. Line types and colors are
as described in Fig. 1. Top row: Different dynamical models for fiducial
values of the C/O ratio, CR ionization rate, and temperature. Seond row:
effect of the CR ionization rate for the fiducial magnetic and non-magnetic
models. Third row: effect of the C/O ratio for each dynamical model, for
fiducial values of the CR ionization rate and temperature. Bottom row: effect
of temperature for each dynamical model for fiducial values of the CR
ionization rate and C/O ratio.  Figure 3.—  As in Fig. 2 but for the HCO+ / NO
pair. Figure 4.—  As in Fig. 2 but for the HCO+ / CO pair. Figure 5.—  As in
Fig. 2 but for the NO+ / NO pair.

Figure 2 shows the results of our aforementioned models for the evolution of
the commonly observed ion-neutral pair HCO+/HCN. We show the evolution of the
central abundance ratio with central density (left column), and radial
profiles of the abundance ratio at two different evolutionary stages (of
central density $10^{4}$ and $10^{6}$ cm-3, central and right columns,
respectively). The top row shows the effect of dynamics on the abundance
ratio, while the second, third, and bottom rows demonstrate the effect of the
CR ionization rate, C/O ratio, and temperature, respectively.

As shown in the left column, the central abundance ratio, especially in the
magnetic models, is not particularly stable with central density. The reason,
as can be seen in Figure 5 of Paper I, is that although both HCO+ and HCN
deplete with increasing central density, in our magnetic models HCO+ depletes
initially faster than HCN and therefore their ratio declines by almost three
orders of magnitude by a central density of $10^{4}$ cm-3. However, at higher
central densities, the depletion rate of HCO+ slows down while HCN keeps
depleting at its initial rate. This effect is responsible for the increase of
the abundance ratio at higher densities. In non-magnetic models, the
dependence of the abundance ratio on central density is much less pronounced,
as the behavior of each of HCO+ and HCN with central density is qualitatively
similar (Figure 5 of Paper I), with only a small rise in the abundance ratio
above $10^{6}$ cm-3 due to some flattening of the depletion rate of HCO+ at
the highest densities examined.

This behavior is also reflected in the radial profiles of this abundance
ratio, as seen in the middle and right columns of Figure 2. The magnetic
models and the slow non-magnetic models (blue, cyan, and purple dashed curves)
exhibit a remarkable variation between the outer layers of the star-forming
core where significant amounts of HCO+ still exist, and the innermost and
densest parts. The fast non-magnetic models, however, have a much milder
dependence on radius. At a more advanced evolutionary stage (n${}_{c}\sim
10^{6}$ cm-3), the minimum seen in the evolution of the central abundance
ratio with density in magnetic models also appears in the radial profile at
intermediate scales. At scales comparable to the size of the core, the radial
profiles of the HCO+/HCN abundance ratio flatten out for all models.

The total variation of the abundance ratio across models, parameters, and
densities is large (several orders of magnitude). As seen in the top row of
Fig. 2, most of the total variation has its origin in the dynamics, with the
CR ionization rate having the second most important contribution (second row).
However, different models represent different possible clouds in nature, and
for this reason the variation of the abundance ratio across models is not the
relevant concern for the Houde et al. type of analysis; instead, the quantity
of interest is the total variation of the ratio within a single model, across
space and time. Still, as seen in for example in the left column of Fig. 2,
the maximum possible variation within a single model can be as large as
several orders of magnitude.

Observationally, the intensity maps of the HCO+/HCN pair show highest degree
of correlation at large length scales, about $0.5-1.0$ pc (Lo et al., 2009),
which correspond to scales comparable with or larger than the entire core in
our models ($r/R_{0}\gtrsim 1$ in Figure 2). However, this correlation can be
reduced at smaller scales of about $0.1$ pc, i.e., resolutions achieved by
interferometric observations (Csengeri et al., 2011).

HCO+ and HCN follow independent chemical networks and thus have no reason to
be co-evolving. The two molecules evolve differently and cannot be considered
co-spatial in most of our models. HCO+ is more abundant at lower densities and
larger spatial scales. Furthermore, the critical densities of the commonly
observed rotational transitions of HCO+ are lower than these of HCN.
Consequently, HCO+ would in principle be expected to trace larger spatial
scales, regions with higher levels of turbulence, which would lead to wider
linewidths for HCO+ rather than HCN. However, this is contradictory to
observations. Houde et al. (2000a) observed the $J=3\rightarrow 2$ and
$J=4\rightarrow 3$ transitions of HCO+ and HCN in several molecular clouds
where they found these species co-existent, with the ion lines indeed more
opaque but narrower than the neutral lines. This can be explained by the
effect of the magnetic field on the dynamics of molecular ions through the
cyclotron interaction, which can in turn be interpreted as a signature of
ambipolar diffusion (Houde et al., 2004) and be used to estimate the strength
of the magnetic field. Indeed, it was later shown that clouds with stronger
magnetic fields (projected in the plane of the sky) show larger differences
between the ion and neutral spectral linewidths (Li et al., 2010; Hezareh et
al., 2012). In contrast, chemical and radiative transfer effects act in the
opposite direction causing the widths of neutral spectral lines to become
narrower than those of ions, thus underestimating the strength of the magnetic
field.

To make a more robust comparison of our models with existing observations we
must perform radiative transfer calculations for the lines observed. We will
return to this problem in the future. An observational correlation similar to
HCO+ and HCN is observed for the H13CO+ and H13CN pair but since our chemical
models do not incorporate the chemistry of the isotopes of carbon and other
elements, we do not include these isotopologue pairs in our discussion. For
the remainder of this paper, we identify additional ion and neutral species
that are good candidates for being considered as co-existent pairs.

###  3.2. HCO+/NO

Figure 3 shows the same pattern as in Figure 2, i.e., the evolution of central
abundance ratio with central density in the left column, and radial profiles
of the abundance ratio for central densities of $10^{4}$ and
$10^{6}{\rm\,cm^{-3}}$ in the middle and right columns, respectively, for the
HCO+ / NO pair. The central abundance ratio of this pair is significantly more
stable with central density within any single model than the HCO+/HCN ratio,
as HCO+ and NO evolve in a similar manner with central density, with the
exception of the high C/O ratio models (third row: yellow and cyan lines for
the magnetic and non-magnetic models, respectively). In all other cases, the
variation of the abundance ratio with central density within a single model is
$\lesssim$ an order of magnitude over the central density range we examine.

A similar behavior is reflected in the radial profiles of the HCO+ / NO
abundance ratio, with most models exhibiting little variation. The large
variation in the HCO+ / NO pair for the magnetic high C/O-ratio model (yellow
curve) only affects very small scales, which are also associated with high
densities. At these densities, however, the depletion of these molecules is
significant (Fig. 5 in Paper I), so these regions are not expected to have a
considerable contribution to the observed lines.

###  3.3. HCO+/CO

The plots for the HCO+ / CO pair is shown in Figure 4. This pair exhibits very
small variations of the central abundance ratio with central density for all
but the fastest non-magnetic models. However, even in the case of these fast
models (dotted blue, cyan, and purple curves), the variation is less than an
order of magnitude for central densities above 104 cm-3.

In the radial profiles, some variation in the abundance ratio occurs at large
radii, but it is less than about an order of magnitude, at least in early
times, and not much more than an order of magnitude in later times, even in
the most-varying, fast non-magnetic models.

It is not surprising that the abundances of the two molecules are co-evolving,
since the dominant reaction for the production of HCO+ in our models involves
CO:

${\rm CO+H_{3}^{+}\rightarrow HCO^{+}+H_{2}}\,,$ (1)

while the destruction of HCO+ proceeds primarily through dissociative
recombination with an electron:

${\rm HCO^{+}+e^{-}\rightarrow CO+H}\,,$ (2)

which leads back to creating CO, forming a production/destruction loop for the
two molecules.

###  3.4. NO+/NO

Figure 5 shows the same plots as before for the NO+ / NO pair. This is the
best example of co-evolution and chemical co-spatiality of an ion-neutral pair
that we have identified in our coupled chemical/dynamical models. The
variation of the central abundance ratio with central density in any single
model is at most a factor of a few. These findings are also reflected in the
radial profiles as well.

The reason for the co-evolution of the two molecules is that NO+ is mainly
produced through ionization of NO or charge-exchange reactions, also involving
NO:

${\rm NO+H^{+}\rightarrow NO^{+}+H}$ (3)

and

${\rm NO+C^{+}\rightarrow NO^{+}+C}\,.$ (4)

The ion is destroyed through dissociative recombination with an electron:

${\rm NO^{+}+e^{-}\rightarrow O+N}\,,$ (5)

so the abundance of NO+ is coupled to that of NO.

##  4\. Discussion

In addition to the co-spatiality of the three new ion-neutral pairs that we
have identified, their observability and critical densities for excitation of
their rotational transitions must also be considered to decide whether they
are viable candidates for studies of the ambipolar ion-neutral drift in star-
forming regions.

The most frequently observed transitions of HCO+ have critical densities111All
molecular data are taken from the Leiden Atomic and Molecular Database
(Schöier et al., 2005), http://home.strw.leidenuniv.nl/$\sim$moldata/  of the
order of 10${}^{5-6}{\rm\,cm^{-3}}$. NO generally has much lower critical
densities, so overall unless higher critical density transitions of NO can be
observed, this would pose a problem in the use of the HCO+/NO pair. The reason
is that NO would trace lower-density, higher-turbulence medium, which would
lead to wider linewidths for NO, a trend that could mimic the effect observed
in the HCO+/HCN pair. Possible higher critical-density transitions would, for
example, be the $J=7/2\rightarrow 5/2$ transition of NO at $350.6$ GHz, with a
critical density of $\simeq 10^{5}{\rm\,cm^{-3}}$. The multiplets of the first
two rotational transitions of NO, i.e., $J=3/2\rightarrow 1/2$ and
$J=5/2\rightarrow 3/2$ have been observed in both dark and massive clouds
(Gerin et al., 1992) at 150 GHz and 250 GHz with relative abundances
calculated to vary from $10^{-7}-10^{-8}$.

A similar problem arises with the HCO+/CO pair, as CO, with generally low
critical densities in usually observed transitions, traces more extended
regions. Because the abundance of CO is high, such transitions would result in
saturated and optically thick lines. To avoid this effect, one needs to use
higher CO transitions. A possibility would be the $J=6\rightarrow 5$
transition at $691$ GHz, with critical density around $10^{5}{\rm\,cm^{-3}}$.

The NO+/NO pair is in principle the best choice in terms of co-spatiality of
the neutral and the ion, but although NO+ has sub-millimeter and radio
transitions identified in the laboratory (Bowman et al., 1982), the ion has
not yet been observed in molecular clouds. This may be due to a general trend
in which nitrogen-bearing molecules, with the exception of NH3 and N2H+, have
not been targeted for observations as much as the carbon-bearing species, in
part because nitrogen chemistry in molecular clouds is even less understood
(especially at low temperatures) than that of carbon (e.g., Hily-Blant et al.,
2010).

We note that our models do not include the effect of turbulence; as a result,
the relative abundance variation that we find is an upper limit, since
turbulent mixing would lessen this variation (Xie et al., 1995). Also, our
models represent prestellar cores, and do not include the effect of feedback
sources, because our parameter studies do not feature very large variations of
the temperature, and UV ionization has not been accounted for. However,
regions with feedback sources are avoided when applying the Houde et al.
method.

We thank Harold Yorke and Helmut Wiesemeyer for insightful and constructive
comments that improved this paper. T.H. is funded by the Alexander von
Humboldt foundation. The project was supported in part by the NASA Origins of
Solar Systems program. Part of this work was carried out at the Jet Propulsion
Laboratory, California Institute of Technology under contract with the
National Aeronautics and Space Administration.

## References

  * Bowman et al. (1982) Bowman, W. C., Herbst, E., & de Lucia, F. C. 1982, J. Chem. Phys., 77, 4261
  * Braginskii (1965) Braginskii, S. I. 1965, Reviews of Plasma Physics, 1, 205
  * Csengeri et al. (2011) Csengeri, T., Bontemps, S., Schneider, N., Motte, F., & Dib, S. 2011, A&A 527, 135 Reviews of Plasma Physics, 1, 205
  * Crutcher et al. (1999) Crutcher, R. M., Troland, T. H., Lazareff, B., Paubert, G., & Kazès, I. 1999, ApJ, 514, 2, L121
  * Gerin et al. (1992) Gerin, M., Viala, Y., Pauzat, F., Ellinger, Y. 1992, A&A, 266, 463
  * Herbst & Klemperer (1973) Herbst, E., & Klemperer, W., 1973, ApJ, 185, 505
  * Hezareh et al. (2012) Hezareh, T., Csengeri, T., Houde, M., Herpin, F., & Bontemps, S. 2012, A&A, submitted
  * Hezareh et al. (2010) Hezareh, T., Houde, M., McCoey, C., & Li, H.-b. 2010, ApJ, 720, 603
  * Hezareh et al. (2008) Hezareh, T., Houde, M., McCoey, C., Vastel, C., & Peng, R 2008, ApJ, 684, 1221
  * Hildebrand et al. (2009) Hildebrand, R. H., Kirby, L., Dotson, J. L., Houde, M., & Vaillancourt, J. E. 2009, ApJ, 696, 567
  * Hily-Blant et al. (2010) Hily-Blant, P., Walmsley, M., Pineau Des Forêts, G., & Flower, D. 2010, A&A, 513, A41
  * Houde et al. (2000a) Houde, M., Bastien, P., Peng, R., Phillips, T. G., & Yoshida, H. 2000a, ApJ, 536, 857
  * Houde et al. (2000b) Houde, M., Peng, R., Phillips, T. G., Bastien, P., & Yoshida, H. 2000, ApJ, 537, 245
  * Houde et al. (2004) Houde, M., Dowell, C. D., Hildebrand, R. H., et al. 2004, ApJ, 604, 717
  * Kulsrud & Pearce (1969) Kulsrud, R., & Pearce, W. P. 1969, ApJ, 156, 445
  * Lai et al. (2003) Lai, S.-P., Velusamy, T., & Langer, W. D. 2003, ApJ, 596, L239
  * Li et al. (2010) Li, H.-b., Houde, M., Lai, S.-p., & Sridharan, T. K. 2010, ApJ, 718, 905
  * Li & Houde (2008) Li, H.-b., & Houde, M. 2008, ApJ, 677, 1151
  * Lo et al. (2009) Lo, N., et al. 2009, MNRAS, 395, 1021
  * McCall et al. (2003) McCall, B. J. 2003, Nature, 422, 500
  * Mestel & Spitzer (1956) Mestel, L., & Spitzer, L., Jr. 1956, MNRAS, 116, 503
  * Mouschovias & Spitzer (1976) Mouschovias, T. C., & Spitzer, L., Jr. 1976, ApJ, 210, 326
  * Mouschovias et al. (2011) Mouschovias, T. C., Ciolek, G. E., & Morton, S. A. 2011, MNRAS, 415, 1751
  * Schöier et al. (2005) Schöier, F. L., van der Tak, F. F. S., van Dishoeck, E. F., & Black, J. H. 2005, A&A, 432, 369
  * Tassis et al. (2012) Tassis, K., Willacy, K., Yorke, H. W., & Turner, N. J. 2012, ApJ, 753, 29, (paper I)
  * Xie et al. (1995) Xie, T., Allen, M., & Langer, W. D. 1995, ApJ, 440, 674
  * Zweibel & Josafatsson (1983) Zweibel, E. G., & Josafatsson, K. 1983, ApJ, 270, 511


# Adversarial Attack and Defense for Non-Parametric Two-Sample Tests

Xilie Xu  Jingfeng Zhang  Feng Liu  Masashi Sugiyama  Mohan Kankanhalli

###### Abstract

Non-parametric two-sample tests (TSTs) that judge whether two sets of samples
are drawn from the same distribution, have been widely used in the analysis of
critical data. People tend to employ TSTs as trusted basic tools and rarely
have any doubt about their reliability. This paper systematically uncovers the
failure mode of non-parametric TSTs through adversarial attacks and then
proposes corresponding defense strategies. First, we theoretically show that
an adversary can upper-bound the distributional shift which guarantees the
attack’s invisibility. Furthermore, we theoretically find that the adversary
can also degrade the lower bound of a TST’s test power, which enables us to
iteratively minimize the test criterion in order to search for adversarial
pairs. To enable TST-agnostic attacks, we propose an ensemble attack (EA)
framework that jointly minimizes the different types of test criteria. Second,
to robustify TSTs, we propose a max-min optimization that iteratively
generates adversarial pairs to train the deep kernels. Extensive experiments
on both simulated and real-world datasets validate the adversarial
vulnerabilities of non-parametric TSTs and the effectiveness of our proposed
defense. Source code is available at https://github.com/GodXuxilie/Robust-
TST.git.

Machine Learning, ICML

##  1 Introduction

Non-parametric two-sample tests (TSTs) that judge whether two sets of samples
drawn from the same distribution have been widely used to analyze critical
data in physics (baldi2014searching), neurophysiology (rasch2008predicting),
biology (borgwardt2006integrating), etc. Compared with traditional methods
(such as the t-test), non-parametric TSTs can relax the strong parametric
assumption about the distributions being studied and are effective in complex
domains (gretton2009fast; gretton2012kernel; chwialkowski2015fast;
jitkrittum2016interpretable; sutherland2016generative; lopez2016revisiting;
cheng2019classification; liu2020learning; liu2021meta). Notably, the use of
deep kernels (liu2020learning) flexibly empowers the non-parametric TSTs to
learn even more complex distributions.

However, the adversarial robustness of non-parametric TSTs is rarely studied,
despite its extensive studies for deep neural networks (DNNs). Studies of
DNNs’ adversarial robustness (Madry_adversarial_training) have enabled
significant advances in defending against adversarial attacks (szegedy), which
can help enhance the security in various domains such as computer vision
(xie2017adversarial; mahmood2021robustness), natural language processing
(Zhu2020FreeLB:; yoo2021towards), recommendation system (peng2020robust), etc.
We therefore undertake this pioneer study on adversarial robustness of non-
parametric TSTs, which uncovers the failure mode of non-parametric TSTs
through adversarial attacks and facilitate an effective strategy for making
TSTs reliable in critical applications (baldi2014searching;
rasch2008predicting; borgwardt2006integrating).

First, we theoretically show the adversary could upper-bound the
distributional shift and degrade the lower bound of a TST’s test power
(details in Section 3.1). Given a benign pair
$(S_{\mathbb{P}},S_{\mathbb{Q}})$, in which
$S_{\mathbb{P}}=\\{{x}_{i}\\}_{i=1}^{m}\sim\mathbb{P}^{m}$ and
$S_{\mathbb{Q}}=\\{{y}_{j}\\}_{j=1}^{n}\sim\mathbb{Q}^{n}$ , an
$\ell_{\infty}$-bounded adversary could generate the adversarial pair
$(S_{\mathbb{P}},\tilde{S}_{\mathbb{Q}})$. We will show in Proposition 1 that
the maximum mean discrepancy (MMD) (gretton2012kernel) between the benign and
adversarial pairs is upper-bounded, which guarantees imperceptible adversarial
perturbations (szegedy). Furthermore, we will show in Theorem
LABEL:theory:tp_attack that the adversary can degrade the lower bound of a
TST’s test power, which implies that a TST could wrongly determine
$\mathbb{P}=\mathbb{Q}$ with a larger probability under adversarial attacks
when $\mathbb{P}\neq\mathbb{Q}$ holds.

Figure 1: An example of adversarial pair
$(S_{\mathbb{P}},\tilde{S}_{\mathbb{Q}})$ generated by embedding an
adversarial perturbation in the benign set $S_{\mathbb{Q}}$ of the benign pair
$(S_{\mathbb{P}},S_{\mathbb{Q}})$. Experimental details are in Section
LABEL:sec:attack_power.

Then, we realize effective adversarial attacks against non-parametric TSTs
(details in Section LABEL:sec:realization). We formulate an attack as a
constraint optimization problem that minimizes a TST’s test criterion
(liu2020learning) within the $\ell_{\infty}$-bound of size $\epsilon$ on
$S_{\mathbb{Q}}$. We utilize projected gradient descent (PGD)
(Madry_adversarial_training) to efficiently search the adversarial set
$\tilde{S}_{\mathbb{Q}}$ and incorporate automatic schedule of the step size
(croce2020reliable) to improve the optimization convergence. Moreover, we
extend the attack beyond a specific TST to a generic TST-agnostic attack,
namely, ensemble attack (EA). EA jointly minimizes a weighted sum of different
test criteria, which can simultaneously fool various TSTs. For example, Figure
1 shows non-parametric TSTs can correctly differentiate the benign pair of
“cats” and “dogs” (top) coming from the different distributions, but wrongly
judge adversarial pairs (bottom) as belonging to the same distribution.

Second, to robustify the non-parametric TSTs, we study the corresponding
defense approaches (details in Section LABEL:sec:defense). A straightforward
defense seems to use an ensemble of TSTs. We find an ensemble of TSTs is
sometimes effective against a specific attack targeting a certain type of TSTs
but almost always fails under EA (see experiments in Section
LABEL:sec:attack_power). Therefore, to effectively defend against adversarial
attacks, we propose to adversarially learn the robust kernels. The defense is
formulated as a max-min optimization that is similar in flavor to the
adversarial training’s min-max formulation (Madry_adversarial_training). For
its realization, we iteratively generate adversarial pairs by minimizing the
test criterion in the inner minimization and update kernel parameters by
maximizing the test criterion on the adversarial pairs in the outer
maximization. We realize our defense using deep kernels that have achieved the
state-of-the-art (SOTA) performance in non-parametric TSTs (liu2020learning).

Lastly, we empirically justify the proposed attacks and defenses (in Section
LABEL:sec:exp). We evaluate the test power of many existing non-parametric
TSTs (non-robust) and the robust-kernel TST (robust) under the EA on simulated
and real-world datasets, including complex synthetic distributions, high-
energy physics data, and challenging images. Comprehensive experimental
results validate that the existing non-parametric TSTs lack adversarial
robustness; we can significantly improve the adversarial robustness of non-
parametric TSTs through adversarially learning the deep kernels.

##  2 Non-Parametric Two-Sample Tests

In this section, we provide the preliminaries of non-parametric TSTs and
provide discussions with the related studies in Appendix
LABEL:sec:related_work.

###  2.1 Problem Formulation

Let $\mathcal{X}\subset\mathbb{R}^{d}$ and $\mathbb{P}$, $\mathbb{Q}$ be Borel
probability measures on $\mathcal{X}$. A non-parametric TST
$\mathcal{J}(S_{\mathbb{P}},S_{\mathbb{Q}}):\mathcal{X}^{m}\times\mathcal{X}^{n}\mapsto\\{0,1\\}$
is used to distinguish between the null hypothesis
$\mathcal{H}_{0}:\mathbb{P}=\mathbb{Q}$ and the alternative hypothesis
$\mathcal{H}_{1}:\mathbb{P}\neq\mathbb{Q}$, where $S_{\mathbb{P}}$ and
$S_{\mathbb{Q}}$ are independent identically distributed (IID) samples of size
$m$ and $n$ drawn from $\mathbb{P}$ and $\mathbb{Q}$, respectively. A non-
parametric TST constructs a mean embedding based on a kernel parameterized
with $\theta$ for each distribution, and utilizes the differences in these
embeddings as the _test statistic_ for the hypothesis test. The judgement is
made by comparing the test statistic
$\mathcal{D}(S_{\mathbb{P}},S_{\mathbb{Q}})$ with a particular threshold $r$:
if the threshold is exceeded, then the test rejects $\mathcal{H}_{0}$. The
_test power_ (TP) of a non-parametric TST $\mathcal{J}$ is measured by the
probability of correctly rejecting $\mathcal{H}_{0}$ when the alternative
hypothesis is true, i.e.,
$\mathrm{TP}(\mathcal{J})=\mathbb{E}_{S_{\mathbb{P}}\sim\mathbb{P}^{m},S_{\mathbb{Q}}\sim\mathbb{Q}^{n}}[\mathbbm{1}(\mathcal{J}(S_{\mathbb{P}},S_{\mathbb{Q}})=1)]$
for a paritular $\mathbb{P}\neq\mathbb{Q}$. A non-parametric TST optimizes its
learnable parameters $\theta$ via maximizing its _test criterion_ , thus
approximately maximizing its test power.

###  2.2 Test Statistics

Here, we introduce a typical test statistic, maximum mean discrepancy (MMD)
(gretton2012kernel), and leave other test statistics in Appendix
LABEL:appendix:tst_intro, such as tests based on Gaussian kernel mean
embeddings at specific positions (chwialkowski2015fast;
jitkrittum2016interpretable) and classifier two-sample tests (C2ST)
(lopez2016revisiting; cheng2019classification).

######  Definition 1 (gretton2012kernel).

Let $k:\mathcal{X}\times\mathcal{X}\to\mathbb{R}$ be a kernel of a reproducing
kernel Hilbert space $\mathcal{H}_{k}$, with feature maps
$k(\cdot,x)\in\mathcal{H}_{k}$. Let $X\sim\mathbb{P}$ and $Y\sim\mathbb{Q}$,
and define the kernel mean embeddings
$\mu_{\mathbb{P}}=\mathbb{E}[k(\cdot,X)]$ and
$\mu_{\mathbb{Q}}=\mathbb{E}[k(\cdot,Y)]$. Under mild integrability
conditions,

$\displaystyle\mathrm{MMD}(\mathbb{P},\mathbb{Q};\mathcal{H}_{k})$
$\displaystyle=\sup_{f\in\mathcal{H},\|f\|_{\mathcal{H}_{k}}\leq
1}|\mathbb{E}[f(X)]-\mathbb{E}[f(Y)]|$
$\displaystyle=\|\mu_{\mathbb{P}}-\mu_{\mathbb{Q}}\|_{\mathcal{H}_{k}}.$ (1)

For characteristic kernels,
$\mathrm{MMD}(\mathbb{P},\mathbb{Q};\mathcal{H}_{k})=0$ if and only if
$\mathbb{P}=\mathbb{Q}$. Assuming $n=m$, we can estimate $\mathrm{MMD}$ (Eq.
(1)) using the $U$-statistic estimator, which is unbiased for
${\mathrm{MMD}}^{2}$ and has nearly minimal variance among all unbiased
estimators (gretton2012kernel):

$\displaystyle{\widehat{\mathrm{MMD}}}^{2}(S_{\mathbb{P}},S_{\mathbb{Q}};k)=\frac{1}{n(n-1)}\sum_{i\neq
j}H_{ij},$ (2) $\displaystyle H_{ij}$
$\displaystyle=k(x_{i},x_{j})+k(y_{i},y_{j})-k(x_{i},y_{j})-k(x_{j},y_{i}),$

where $x_{i},x_{j}\in S_{\mathbb{P}}$ and $y_{i},y_{j}\in S_{\mathbb{Q}}$.

In this paper, we investigate six types of non-parametric TSTs as follows
since liu2020learning; liu2021meta have shown they are powerful on complex
data.

  * •

$\mathcal{D}^{(\mathrm{G})}(\cdot,\cdot;k^{(\mathrm{G})})={\widehat{\mathrm{MMD}}}^{2}(\cdot,\cdot;k^{(\mathrm{G})})$
for tests based on MMD with Gaussian kernels (MMD-G)
(sutherland2016generative) with the learnable lengthscale $\sigma_{\phi}$, in
which $k^{(\mathrm{G})}(x,y)=\exp(-\frac{1}{2\sigma_{\phi}}\|x-y\|^{2})$.

  * •

$\mathcal{D}^{(\mathrm{D})}(\cdot,\cdot;k^{(\mathrm{D})})={\widehat{\mathrm{MMD}}}^{2}(\cdot,\cdot;k^{(\mathrm{D})})$
for tests based on MMD with deep kernels (MMD-D) (liu2020learning). Note that
$k^{(\mathrm{D})}(x,y)=[(1-\gamma)\exp(-\frac{1}{2\sigma_{\phi}}\|\phi(x)-\phi(y)\|^{2})+\gamma]\exp(-\frac{1}{2\sigma_{q}}\|x-y\|^{2})$
where $\gamma,\sigma_{\phi},\sigma_{q}$ are the learnable parameters and
$\phi(\cdot)$ is a parameterized deep network to extract the features.

  * •

$\mathcal{D}^{(\mathrm{S})}(\cdot,\cdot)$ (Eq. (LABEL:eq:ts_c2sts)) for C2ST
based on Sign (C2ST-S) (lopez2016revisiting). A classifier
$f:\mathcal{X}\rightarrow\mathbb{R}$ that outputs the classification
probabilities is utilized by C2ST. liu2020learning pointed out that the test
statistic of C2ST-S is equivalent to MMD with kernel $k^{(\mathrm{S})}$, i.e.,
$\mathcal{D}^{(\mathrm{S})}(\cdot,\cdot)={\widehat{\mathrm{MMD}}}^{2}(\cdot,\cdot;k^{(\mathrm{S})})$
where
$k^{(\mathrm{S})}(x,y)=\frac{1}{4}\mathbbm{1}(f(x)>0)\mathbbm{1}(f(y)>0)$.

  * •

$\mathcal{D}^{(\mathrm{L})}(\cdot,\cdot)$ (Eq. (LABEL:eq:ts_s2stl)) for C2ST-L
(cheng2019classification) that utilizes the discriminator’s measure of
confidence. Its test statistic is also equivalent to MMD with kernel
$k^{(\mathrm{L})}$ (liu2020learning), i.e.,
$\mathcal{D}^{(\mathrm{L})}(\cdot,\cdot)={\widehat{\mathrm{MMD}}}^{2}(\cdot,\cdot;k^{(\mathrm{L})})$
where $k^{(\mathrm{L})}(x,y)=f(x)f(y)$.

  * •

$D^{(\mathrm{ME})}(\cdot,\cdot)$ (Eq. (LABEL:eq:ts_ME)) for tests based on
differences in Gaussian kernel mean embeddings at specific locations
(chwialkowski2015fast; jitkrittum2016interpretable), namely Mean Embedding
(ME).

  * •

$D^{(\mathrm{SCF})}(\cdot,\cdot)$ (Eq. (LABEL:eq:ts_SCF)) for tests based on
Gaussian kernel mean embeddings at a set of optimized frequency
(chwialkowski2015fast; jitkrittum2016interpretable), namely Smooth
Characteristic Functions (SCF).

###  2.3 Test Criterion

In this subsection, we introduce the test criteria for non-parametric TSTs
based on MMD (sutherland2016generative; liu2020learning; lopez2016revisiting;
cheng2019classification).

######  Theorem 1 (Asymptotics of MMD under $\mathcal{H}_{1}$
(serfling2009approximation)).

Under the alternative, $\mathcal{H}_{1}:\mathbb{P}\neq\mathbb{Q}$, the
standard central limit theorem holds:

$\displaystyle\sqrt{n}({\widehat{\mathrm{MMD}}}^{2}-{\mathrm{MMD}}^{2})\rightarrow\mathcal{N}(0,\sigma_{\mathcal{H}_{1}}^{2}),$
$\displaystyle\sigma_{\mathcal{H}_{1}}^{2}=4(\mathbb{E}[H_{12}H_{13}]-\mathbb{E}[H_{12}]^{2}),$

where $H_{12},H_{13}$ refer to $H_{ij}$ in Eq. (2).

Guided by the asymptotics of MMD (Theorem 1), the test power is estimated as
follows:

$\displaystyle\Pr(n{\widehat{\mathrm{MMD}}}^{2}>r)\rightarrow\Phi(\frac{\sqrt{n}{\mathrm{MMD}}^{2}}{\sigma_{\mathcal{H}_{1}}}-\frac{r}{\sqrt{n}\sigma_{\mathcal{H}_{1}}}),$
(3)

where $\Phi$ is the cumulative distribution function (CDF) of standard normal
distribution and $r$ is the rejection threshold approximately found via
permutation testing (dwass1957modified; fernandez2008test). This general
method is usually considered best to estimate the null hypothesis: under
$\mathcal{H}_{0}$, samples from $\mathbb{P}$ and $\mathbb{Q}$ are
interchangeable, and repeatedly re-computing the test statistic with samples
randomly shuffled between $S_{\mathbb{P}}$ and $S_{\mathbb{Q}}$ estimates its
null distribution.

For reasonably large $n$, the test power is dominated by the first term of Eq.
(3), and thus the TST yields the most powerful test by approximately
maximizing the test criterion (liu2020learning)

$\displaystyle\mathcal{F}(\mathbb{P},\mathbb{Q};k)={\mathrm{MMD}}^{2}(\mathbb{P},\mathbb{Q};k)/\sigma_{\mathcal{H}_{1}}(\mathbb{P},\mathbb{Q};k).$
(4)

Further, $\mathcal{F}(\mathbb{P},\mathbb{Q};k)$ can be empirically estimated
with

$\displaystyle\hat{\mathcal{F}}(S_{\mathbb{P}},S_{\mathbb{Q}};k)=\frac{{\widehat{\mathrm{MMD}}}^{2}(S_{\mathbb{P}},S_{\mathbb{Q}};k)}{\hat{\sigma}_{\mathcal{H}_{1},\lambda}(S_{\mathbb{P}},S_{\mathbb{Q}};k)},$
(5)

where $\hat{\sigma}_{\mathcal{H}_{1},\lambda}^{2}$ is a regularized estimator
of $\sigma_{\mathcal{H}_{1}}^{2}$:

$\displaystyle\hat{\sigma}_{\mathcal{H}_{1},\lambda}^{2}=\frac{4}{n^{3}}\sum_{i=1}^{n}\bigg{(}\sum_{j=1}^{n}H_{ij}\bigg{)}^{2}-\frac{4}{n^{4}}\bigg{(}\sum_{i=1}^{n}\sum_{j=1}^{n}H_{ij}\bigg{)}^{2}+\lambda,$

where $\lambda$ is a positive constant. The test criterion of the MMD test
(e.g., MMD-G, MMD-D, C2ST-S and C2ST-L) is calculated based on its
corresponding kernel. We let $\sigma_{\theta}^{2}$ denote
$\sigma^{2}_{\mathcal{H}_{1}}(\mathbb{P},\mathbb{Q};k_{\theta})$, and
analogously $\hat{\sigma}_{\theta}^{2}$ denote
$\hat{\sigma}^{2}_{\mathcal{H}_{1},\lambda}(\mathbb{P},\mathbb{Q};k_{\theta})$,
for simplicity.

In addition, chwialkowski2015fast and jitkrittum2016interpretable analyzed
that the test power of ME tests, and SCF tests can be approximately maximized
by maximizing the corresponding test criterion as well, i.e.,
$\hat{\mathcal{F}}^{(\mathrm{ME})}(S_{\mathbb{P}},S_{\mathbb{Q}})$ and
$\hat{\mathcal{F}}^{(\mathrm{SCF})}(S_{\mathbb{P}},S_{\mathbb{Q}})$ (details
in Appendix LABEL:appendix:tst_intro).

To avoid notation clutter, we simply let $\theta$ represent all of the
learnable parameters in a non-parametric TST. The optimized parameters of a
TST are obtained as follows.

$\displaystyle\hat{\theta}\approx\operatorname*{\mathrm{arg\,max}}_{\theta}\hat{\mathcal{F}}(S_{\mathbb{P}}^{\mathrm{tr}},S_{\mathbb{Q}}^{\mathrm{tr}};k_{\theta}),$
(6)

where $(S_{\mathbb{P}}^{\mathrm{tr}},S_{\mathbb{Q}}^{\mathrm{tr}})$ is the
training pair. Then, we conduct a hypothesis test based on
$\mathcal{D}(S_{\mathbb{P}}^{\mathrm{te}},S_{\mathbb{Q}}^{\mathrm{te}};k_{\hat{\theta}})$,
where $(S_{\mathbb{P}}^{\mathrm{te}},S_{\mathbb{Q}}^{\mathrm{te}})$ is the
test pair.

##  3 Adversarial Attacks Against Non-Parametric TSTs

In this section, we first show the possible existence of adversarial attacks
against a non-parametric TST. Then, we propose a method to generate
adversarial test pairs that can fool a TST. To enable TST-agnostic attacks, we
propose a unified attack framework, i.e., ensemble attack.

###  3.1 Theoretical Analysis

This section theoretically shows that there could exist adversarial attacks
that can invisibly undermine a TST. We first lay out the needed assumptions on
kernel functions.

######  Assumption 1.

The possible kernel parameterized with $\theta\in\mathbb{R}^{\kappa}$ lies in
Banach space. The set of possible kernel parameters $\Theta$ is bounded by
$\mathrm{R}_{\Theta}$, i.e.,
$\Theta\subseteq\\{\theta\mid\|\theta\|\leq\mathrm{R}_{\Theta}\\}$. We let
$\bar{\Theta}_{s}=\\{\theta\in\Theta\mid\sigma_{\theta}^{2}\geq s^{2}>0\\}$ in
which $s$ is a positive constant.

######  Assumption 2.

The kernel function $k_{\theta}$ is uniformly bounded, i.e.,
$\sup_{\theta\in\Theta}\sup_{x\in\mathcal{X}}k_{\theta}(x,x)\leq\nu$. We treat
$\nu$ as a constant.

######  Assumption 3.

The kernel function $k_{\theta}(x,y)$ satisfies the Lipschitz conditions as
follows.

$\displaystyle|k_{\theta}(x,y)-k_{\theta^{\prime}}(x,y)|\leq
L_{1}\|\theta-\theta^{\prime}\|;$
$\displaystyle|k_{\theta}(x,y)-k_{\theta}(x^{\prime},y^{\prime})|\leq
L_{2}(\|x-x^{\prime}\|+\|y-y^{\prime}\|),$

where $L_{1}$ and $L_{2}$ are positive constants.

We consider a potential risk that causes a malfunction of a non-parametric
TST: an adversarial attacker that aims to deteriorate the TST’s test power,
can craft an adversarial pair $(S_{\mathbb{P}},\tilde{S}_{\mathbb{Q}})$ as the
input to the TST during the testing procedure, in which the two sets
$\tilde{S}_{\mathbb{Q}}$ and $S_{\mathbb{Q}}$ are nearly indistinguishable. We
provide a detailed description of the attacker against non-parametric TSTs in
Appendix LABEL:appendix:attacker.

We define the $\epsilon$-ball centered at $x\in\mathcal{X}$ as follows:

$\displaystyle\mathcal{B}_{\epsilon}[x]=\\{\tilde{x}\in\mathcal{X}\mid\|x-\tilde{x}\|_{\infty}\leq\epsilon\\}.$

Further, an $\ell_{\infty}$-bound of size $\epsilon$ on the set
$S_{\mathbb{Q}}$ is defined as

$\displaystyle\mathcal{B}_{\epsilon}[S_{\mathbb{Q}}]=\\{$
$\displaystyle\tilde{S}_{\mathbb{Q}}=\\{\tilde{x}_{i}\in\mathcal{X}\\}_{i=1}^{n}\mid$
$\displaystyle\tilde{x}_{i}\in\mathcal{B}_{\epsilon}[x_{i}],\quad\forall
x_{i}\in S_{\mathbb{Q}},\tilde{x}_{i}\in\tilde{S}_{\mathbb{Q}}\\}.$

Without loss of generality, we assume that the adversarial perturbation is
$\ell_{\infty}$-bounded of size $\epsilon$, i.e.,
$\tilde{S}_{\mathbb{Q}}\in\mathcal{B}_{\epsilon}[S_{\mathbb{Q}}]$. We leave
exploring the effects of other constraints that can bound the “human
imperception” as the future work, such as Wasserstein-distance constraints
(wong2019wasserstein).

Under $\ell_{\infty}$-bounded attacks, we conduct our theoretical analysis of
distributional shift in the test pairs as follows.

######  Proposition 1.

Under Assumptions 1 to 3, we use $n_{\mathrm{tr}}$ samples to train a kernel
$k_{\theta}$ parameterized with $\theta$ and $n_{\mathrm{te}}$ samples to run
a test of significance level $\alpha$. Given the adversarial budget
$\epsilon\geq 0$, the benign pair $(S_{\mathbb{P}},S_{\mathbb{Q}})$ and the
corresponding adversarial pair $(S_{\mathbb{P}},\tilde{S}_{\mathbb{Q}})$ where
$\tilde{S}_{\mathbb{Q}}\in\mathcal{B}_{\epsilon}[S_{\mathbb{Q}}]$, with the
probability at least $1-\delta$, we have

$\displaystyle\sup_{\theta}|{\widehat{\mathrm{MMD}}}^{2}(S_{\mathbb{P}},\tilde{S}_{\mathbb{Q}};k_{\theta})-{\widehat{\mathrm{MMD}}}^{2}(S_{\mathbb{P}},S_{\mathbb{Q}};k_{\theta})|$
$\displaystyle\leq\frac{8L_{2}\epsilon\sqrt{d}}{\sqrt{n_{\mathrm{te}}}}\sqrt{2\log\frac{2}{\delta}+2\kappa\log(4\mathrm{R}_{\Theta}\sqrt{n_{\mathrm{te}}})}+\frac{8L_{1}}{\sqrt{n_{\mathrm{te}}}}.$

The proof is in Appendix LABEL:appendix:proof_mmd_shift.

######  Remark 1.

Proposition 1 shows that $\epsilon$ can control the upper bound of
distributional shift measured by MMD between samples in the test pair. In
other words, a small $\epsilon$ can ensure the difference between
${\widehat{\mathrm{MMD}}}^{2}(S_{\mathbb{P}},\tilde{S}_{\mathbb{Q}};k_{\theta})$
and ${\widehat{\mathrm{MMD}}}^{2}(S_{\mathbb{P}},S_{\mathbb{Q}};k_{\theta})$
is numerically small. Therefore, an $\ell_{\infty}$-bounded adversary can make
the adversarial perturbation imperceptible, thus guaranteeing the attack’s
invisibility.

Next, we provide a lemma that theoretically analyzes the adversary’s influence
on the estimated test criterion.

######  Lemma 1.

In the setup of Proposition 1, with probability at least $1-\delta$, we have

$\displaystyle\sup_{\theta\in\bar{\Theta}_{s}}|\hat{\mathcal{F}}(S_{\mathbb{P}},\tilde{S}_{\mathbb{Q}};k_{\theta})-\hat{\mathcal{F}}(S_{\mathbb{P}},S_{\mathbb{Q}};k_{\theta})|$
$\displaystyle=\mathcal{O}\bigg{(}\frac{\epsilon
L_{2}\sqrt{d\big{(}\log\frac{1}{\delta}+\kappa\log(\mathrm{R}_{\Theta}\sqrt{n_{\mathrm{te}}})\big{)}}+L_{1}}{s\sqrt{n_{\mathrm{te}}}}\bigg{)}.$

The proof is in Appendix LABEL:appendix:proof_adv_benign_tp.