# SCAT Uncovers ATLAS’s First Tidal Disruption Event ATLAS18mlw: A Faint and
Fast TDE in a Quiescent Balmer Strong Galaxy

Jason T. Hinkle Institute for Astronomy, University of Hawai‘i, 2680 Woodlawn
Drive, Honolulu, HI 96822, USA  Michael A. Tucker DOE CSGF Fellow  Institute
for Astronomy, University of Hawai‘i, 2680 Woodlawn Drive, Honolulu, HI 96822,
USA  Benjamin. J. Shappee Institute for Astronomy, University of Hawai‘i, 2680
Woodlawn Drive, Honolulu, HI 96822, USA  Thomas W.-S. Holoien NHFP Einstein
Fellow  The Observatories of the Carnegie Institution for Science, 813 Santa
Barbara Street, Pasadena, CA 91101, USA  Patrick J. Vallely Department of
Astronomy, The Ohio State University, 140 West 18th Avenue, Columbus, OH
43210, USA  Thomas de Jaeger  Institute for Astronomy, University of Hawai‘i,
2680 Woodlawn Drive, Honolulu, HI 96822, USA  Katie Auchettl School of
Physics, The University of Melbourne, Parkville, VIC 3010, Australia  ARC
Centre of Excellence for All Sky Astrophysics in 3 Dimensions (ASTRO 3D)
Department of Astronomy and Astrophysics, University of California, Santa
Cruz, CA 95064, USA  Greg Aldering  Lawrence Berkeley National Laboratory, 1
Cyclotron Rd., Berkeley, CA, 94720, USA  Chris Ashall  Institute for
Astronomy, University of Hawai‘i, 2680 Woodlawn Drive, Honolulu, HI 96822, USA
Dhvanil D. Desai Institute for Astronomy, University of Hawai‘i, 2680 Woodlawn
Drive, Honolulu, HI 96822, USA  Aaron Do  Institute for Astronomy, University
of Hawai‘i, 2680 Woodlawn Drive, Honolulu, HI 96822, USA  Anna V. Payne NASA
Fellowship Activity Fellow  Institute for Astronomy, University of Hawai‘i,
2680 Woodlawn Drive, Honolulu, HI 96822, USA  John L. Tonry  Institute for
Astronomy, University of Hawai‘i, 2680 Woodlawn Drive, Honolulu, HI 96822, USA

###### Abstract

We present the discovery that ATLAS18mlw was a tidal disruption event (TDE) in
the galaxy WISEA J073544.83+663717.3, at a distance of 334 Mpc. Initially
discovered by the Asteroid Terrestrial Impact Last Alert System (ATLAS) on
2018 March 17.3, the TDE nature of the transient was uncovered only recently
with the re-reduction of a SuperNova Integral Field Spectrograph (SNIFS)
spectrum. This spectrum, taken by the Spectral Classification of Astronomical
Transients (SCAT) survey, shows a strong blue continuum and a broad H$\alpha$
emission line. Here we present roughly six years of optical survey photometry
beginning before the TDE to constrain AGN activity, optical spectroscopy of
the transient, and a detailed study of the host galaxy properties through
analysis of archival photometry and a host spectrum. ATLAS18mlw was detected
in ground-based light curves for roughly two months. From a blackbody fit to
the transient spectrum and bolometric correction of the optical light curve,
we conclude that ATLAS18mlw is likely a low-luminosity TDE with a peak
luminosity of log(L [erg s-1]) = $43.5\pm 0.2$. The TDE classification is
further supported by the quiescent Balmer strong nature of the host galaxy. We
also calculated the TDE decline rate from the bolometric light curve and find
$\Delta L_{40}=-0.7\pm 0.2$ dex, making ATLAS18mlw a member of the growing
class of “faint and fast” TDEs with low peak luminosities and fast decline
rates.

Black hole physics (159) — Galaxy properties(615) — Supermassive black holes
(1663) — Tidal disruption (1696) — Transient sources (1851)

††journal: ApJ††facilities: ASAS-SN (Shappee et al., 2014; Kochanek et al.,
2017), ATLAS (Tonry et al., 2018), GMOS (Hook et al., 2004), SNIFS (Lantz et
al., 2004), ZTF (Bellm et al., 2019)††software: Linmix (Kelly, 2007),
Matplotlib (Hunter, 2007), Numpy (Harris et al., 2020), PYPHOT, pPXF
(Cappellari & Emsellem, 2004; Cappellari, 2017), Scipy (Virtanen et al., 2020)

##  1 Introduction

When a star passes within the tidal radius of a supermassive black hole
(SMBH), the self-gravity of the star is overwhelmed by tidal forces and the
star is disrupted, resulting in a tidal disruption event (TDE; Rees, 1988;
Evans & Kochanek, 1989; Phinney, 1989). A TDE results in a luminous flare as
the disrupted stellar material falls back onto the SMBH (e.g., Ulmer, 1999).
To date, $\sim 60$ of these flares have been observed (Gezari, 2021), with a
majority discovered in the optical. This is in large part due to the recent
expansion of transient surveys such as the All-Sky Automated Survey for
Supernovae (ASAS-SN; Shappee et al., 2014; Kochanek et al., 2017), the
Asteroid Terrestrial Impact Last Alert System (ATLAS; Tonry et al., 2018),
Gaia Alerts (Wyrzykowski et al., 2012), the Panoramic Survey Telescope and
Rapid Response System (Pan-STARRS; Chambers et al., 2016), and the Zwicky
Transient Facility (ZTF; Bellm et al., 2019).

However, it is not only optical sky surveys that play an important role in
discovering TDEs. Spectroscopic classification surveys are crucial in
obtaining rapid spectra to not only classify transients but to also identify
interesting sources. Alongside the expansion of transient surveys has been the
growth of such spectroscopic surveys. In the past two years, roughly three-
quarters of all spectroscopically-confirmed transients have been classified by
the Zwicky Transient Facility classification efforts (e.g., Fremling et al.,
2020, 2021), the Public ESO Spectroscopic Survey for Transient Objects
(PESSTO; Smartt et al., 2015) and successors (e.g., Barbarino et al., 2019),
and the Spectral Classification of Astronomical Transients (SCAT) survey
(Tucker et al., in prep.), in that order.

Spectroscopic classification of TDEs has become easier as our understanding of
the typical observed emission from TDEs becomes more complete. The emission of
TDEs is hot, typically exhibiting a blue optical continuum, with a UV excess.
The UV/optical spectral energy distribution (SED) of TDEs is often well-fit by
a blackbody with a temperature on the order of T $=20000-50000$ K (e.g.,
Holoien et al., 2014a, 2016; Blagorodnova et al., 2017; Hinkle et al., 2021a).
Roughly half of optically-selected TDEs show X-ray emission, commonly in the
form of a blackbody with kT $\sim 50$ eV, or a soft power-law with a photon
index of $\gtrsim 2.5$ (Brown et al., 2017; Auchettl et al., 2017, 2018). The
characteristic lifetime of a TDE in the optical bands is on the order of a few
months (e.g., Hinkle et al., 2020; van Velzen et al., 2020) although UV
emission can exist above host galaxy levels for several years after disruption
(e.g., van Velzen et al., 2019).

TDEs represent an excellent opportunity to study otherwise non-active SMBHs.
For example, TDE emission is likely sensitive to black hole spin and mass
(e.g. Ulmer, 1999; Graham et al., 2001; Mockler et al., 2019; Gafton &
Rosswog, 2019). It has been shown that the SMBH masses obtained through fits
to TDE light curves are consistent with those derived through other means
(Mockler et al., 2019). This is particularly important, as the large majority
of SMBHs are not active ($\sim 90$%; Ho, 2008; Lacerda et al., 2020).

Unfortunately, TDEs are rare, with observational rates estimated to be between
$10^{-4}$ and $10^{-5}\text{ yr}^{-1}$ per galaxy (e.g. van Velzen & Farrar,
2014; Holoien et al., 2016; van Velzen, 2018; Auchettl et al., 2018).
Interestingly, TDEs seem to prefer particular host galaxies where a recent
burst of star formation has occurred (e.g., Stone & van Velzen, 2016; French
et al., 2016) and moderately massive stars may be preferred (Mockler et al.,
2022). These include quiescent Balmer strong (QBS) and post-starburst (PSB)
galaxies, both likely post-merger systems (French, 2021). TDEs occur at rates
enhanced by factors of roughly 20 and 40 times the typical values for QBS and
PSB galaxies respectively (e.g. Arcavi et al., 2014; French et al., 2016; Law-
Smith et al., 2017; Graur et al., 2018; French et al., 2020b).

Due in part to the preference for TDEs to occur in unique galaxies and their
slightly increased rates around lower mass SMBHs (Wang & Merritt, 2004; Stone
& Metzger, 2016; van Velzen, 2018), the host galaxies of TDEs have received
special attention (Law-Smith et al., 2017; French et al., 2020a, b). To
fulfill the promise of TDEs as probes of SMBHs, studies of the history of
accretion, mergers, and star-formation in the TDE host galaxies are needed
(e.g., Prieto et al., 2016).

Here we present the observations of the first TDE discovered by ATLAS,
ATLAS18mlw. The paper is organized in the following manner. In Section 2 we
detail observations of the TDE candidate and host galaxy. In Section 3 we
discuss the properties of the host galaxy of ATLAS18mlw. In Section 4 we
analyze the data on this source, and in Section 5 we discuss our results.
Throughout this paper, we have used a cosmology with $H_{0}$ = 69.6 km s-1
Mpc-1, $\Omega_{M}=0.29$, and $\Omega_{\Lambda}=0.71$ (Wright, 2006; Bennett
et al., 2014).

##  2 Data

###  2.1 Archival Photometry

We obtained archival photometry of the host galaxy WISEA J073544.83+663717.3
including $ugriz$ magnitudes from the Sloan Digital Sky Survey (SDSS) Data
Release 15 (Aguado et al., 2019) and $W1$ and $W2$ magnitudes from the Wide-
field Infrared Survey Explorer (WISE; Wright et al., 2010) AllWISE catalog.
While the host galaxy was not listed in any Galaxy Evolution Explorer (GALEX;
Martin et al., 2005) catalog, we measured a 6$\farcs$0 radius aperture $FUV$
and $NUV$ limit and magnitude from GALEX images using gPhoton (Million et al.,
2016). These ultraviolet, optical, and infrared host-galaxy magnitudes and
estimated uncertainties are presented in Table 1.

Table 1: Archival Host-Galaxy Photometry

Filter | Magnitude | Magnitude Uncertainty
---|---|---
$FUV$ |  $>$23.56 | 3$\sigma$
$NUV$ | 23.66 | 0.34
$u$ | 21.76 | 0.23
$g$ | 20.08 | 0.03
$r$ | 19.34 | 0.02
$i$ | 18.96 | 0.02
$z$ | 18.79 | 0.06
$W1$ | 19.04 | 0.06
$W2$ | 19.46 | 0.16

Archival magnitudes of the host galaxy WISEA J073544.83+663717.3. The $FUV$
and $NUV$ magnitudes were measured from GALEX images, the $ugriz$ magnitudes
are taken from the SDSS catalog and the $W1$ and $W2$ magnitudes are taken
from the AllWISE catalog. All magnitudes are presented in the AB system.

###  2.2 Transient Imaging

Figure 1: Host-subtracted and Galactic extinction-corrected light curve of the
host galaxy WISEA J073544.83+663717.3 from ATLAS, ZTF, and TESS spanning
approximately six years. The ATLAS18mlw flare is clearly visible during early
2018. The green and red circles represent ZTF $g$\- and $r$-band photometry
stacked in 10-day bins. The cyan and orange circles represent ATLAS $c$\- and
$o$-band photometry stacked in 10-day bins. The dark red circles represent
TESS photometry stacked in 1-day bins. The vertical gray line marks the epoch
of ATLAS discovery.

The transient ATLAS18mlw, $(\alpha,\delta)=$ (07:35:44.777, $+$66:37:16.43),
was discovered on 2018 March 17.3 by the ATLAS survey. The discovery was
announced publicly on the Transient Name Server (TNS) and given the name
AT2018ahl111https://www.wis-tns.org/object/2018ahl. However, rather than
obscure the discovering survey, we will refer to the transient by its survey
name ATLAS18mlw throughout the remainder of this manuscript.

Soon after the ATLAS discovery, we obtained a $V$-band acquisition image with
the Supernova Integral Field Spectrograph (SNIFS; Lantz et al., 2004) on the
University of Hawai‘i 88-in telescope (UH88). From this, we are able to
compare the location of the transient to the location of the host galaxy to
determine if the source is spatially consistent with the nucleus. We first
used astrometry.net (Barron et al., 2008; Lang et al., 2010) to compute
astrometry for both the SNIFS acquistion image and a Pan-STARRS $r$-band image
(Chambers et al., 2016). We then measured the centroid of the host galaxy from
the Pan-STARRS image and the centroid of the TDE in the SNIFS image with the
IRAF task imcentroid. We find a transient position of
$(\alpha,\delta)=($07:35:44.793$,+$66:37:16.53$)$ and a host galaxy position
of $(\alpha,\delta)=($07:35:44.769$,+$66:37:16.34$)$. Combined, these give an
offset of 0$\farcs$4$\pm$0$\farcs$5, where the uncertainty is dominated by the
uncertainties in the astrometric solutions for the SNIFS and Pan-STARRS
images. At the host distance, this offset corresponds to a physical distance
of $700\pm 900$ pc, consistent with ATLAS18mlw being a nuclear transient.

We then obtained the ATLAS light curves of ATLAS18mlw from the ATLAS forced
point-spread function (PSF) photometry service. ATLAS uses two 0.5-m
telescopes located on Haleakalā and at the Mauna Loa Observatory. During
typical operation, the ATLAS telescopes obtain four 30 s exposures covering
roughly a quarter of the sky per night (Smith et al., 2020). ATLAS obtains
images in two broad-band filters, the ‘cyan’ ($c$) filter between 420–650 nm
and the ‘orange’ ($o$) filter from 560–820 nm (Tonry et al., 2018). We combine
intra-night epochs using a weighted mean to obtain deeper limits and more
robust detections.

In addition to ATLAS photometry, we obtained ZTF photometry in the $g$\- and
$r$-bands from their forced PSF photometry service. The ZTF survey makes use
of the Samuel Oschin 48-in telescope at the Palomar Observatory. In normal
operation, ZTF images down to a depth of $\sim$20.5 $r$-band mag in a 30 s
exposure (Bellm et al., 2019). Similar to the ATLAS data, we combined the
intra-night photometric observations using a weighted average to get a single
flux measurement.

We also obtained a light curve of the host galaxy from the Transiting
Exoplanet Survey Satellite (TESS; Ricker et al., 2015). TESS observed the
location of ATLAS18mlw during Sector 20, at a phase of approximately 650 days
after discovery. The TESS full frame images (FFIs) were reduced using the ISIS
package (Alard & Lupton, 1998; Alard, 2000). We performed image subtraction
following the procedures of (Vallely et al., 2019, 2021) to compute our
differential light curves. The reference image was constructed from the first
100 FFIs of good quality. The measured count rates were converted into fluxes
using the instrumental zero point of 20.44 electrons per second (Vanderspek et
al., 2018). To quantify the host variability, we binned the TESS data in 12
hour bins and computed the root mean squared scatter, finding a value of 8.1
$\mu$Jy, which indicates no significant host-galaxy variability in quiescence.

Finally, we obtained ASAS-SN photometry in the $g$\- and $V$-bands from the
newly-public subtracted light curve server222https://asas-sn.osu.edu/
(Kochanek et al., 2017). Over the full baseline of the ASAS-SN survey, the
host galaxy of ATLAS18mlw was imaged 285 times in $V$ and 349 in $g$ with no
images showing significant detections. The median 5$\sigma$ upper limits are
17.6 and 18.0 mag for the $V$\- and $g$-bands respectively.

Figure 1 shows the ATLAS, ZTF, and TESS data for ATLAS18mlw both during the
flare and in periods of quiescence. Here the ATLAS and ZTF data have been
stacked in 10 day bins and the TESS data has been stacked in 1 day bins.
Besides the ATLAS18mlw flare itself, it is apparent that no other large flare
exists in the roughly six years of coverage. This suggests that the host
galaxy did not exhibit AGN-like activity prior to the detection of ATLAS18mlw.

The subset of the light curve during the ATLAS18mlw flare is shown in Figure
2, showing that ATLAS18mlw rose to peak within $\sim 10$ days after discovery,
although there are no non-detections near the discovery epoch to help
constrain the time of first light further. After peak, the transient fades by
two mags within sixty days and is not detected afterwards. This is faster than
most TDEs, although there is no simultaneous UV coverage to test the longevity
of the UV excess commonly seen in TDEs (e.g., Holoien et al., 2016; van Velzen
et al., 2020; Gezari, 2021).

Figure 2: Host-subtracted and Galactic extinction-corrected optical light
curves of ATLAS18mlw. Shown are ATLAS ($co$) and ZTF ($gr$) photometry stacked
in 1-day bins. Downward-facing triangles indicate 3$\sigma$ upper limits in
filters of the same respective color. The thick black bar on the bottom axis
marks the epoch of the initial SNIFS spectrum. All data are shown in the AB
magnitude system.

###  2.3 Spectroscopy

A classification spectrum was obtained 1.05 days after discovery by the SCAT
survey (Tucker et al., 2018) using SNIFS on UH88. SNIFS is comprised of a blue
(B; 320 – 560 nm) and red (R; 520 nm – 1000 nm) channel. The original
classification spectrum was reduced with a summit-based quick-look pipeline
which is a simplified version of the full SNfactory pipeline (Aldering et al.,
2006). The quick-look pipeline estimates the extraction location by averaging
over the wavelength axis in each channel and then does aperture photometry.

The quick-look pipeline is reliable when the wavelength dependence of the
trace position is less than the seeing. This is true for targets observed near
zenith, but atmospheric dispersion shifts the object’s location on the chip as
the airmass increases. ATLAS18mlw was observed at an airmass of $\approx 1.62$
which introduced a wavelength-dependence to the relative contributions from
the TDE and the host-galaxy to the final spectrum, masking the transient’s
signal and hindering classification.

The development of an upgraded reduction pipeline at the University of Hawai‘i
in June 2021 incorporates wavelength-dependent tracing across each channel.
Reanalysis of the classification spectrum revealed a strong blue continuum
with narrow host lines and a broad H$\alpha$ line at a redshift of z = 0.073.
These spectral features, along with the peak absolute magnitude of $-18.6$
mag, are consistent with a tidal disruption event (e.g., van Velzen et al.,
2020). A full description of the updated reduction pipeline will be presented
in a forthcoming publication (Tucker et al., in prep.).

After concluding that ATLAS18mlw was likely a TDE, we obtained a spectrum of
the host galaxy WISEA J073544.83+663717.3 with the Gemini Multi-Object
Spectrograph (GMOS; Hook et al., 2004) on the 8.1-m Gemini North telescope.
This spectrum was taken with the B600 grating and dithered both along the
spatial and spectral directions to mitigate the presence of chip gaps and
avoid the influence of bad pixels. The spectra were reduced using normal
methods such as bias, flat-fielding, and wavelength calibration using arc
lamps. The individual dithered spectra were median combined to minimize the
effect of cosmic ray hits. Figure 3 shows the SNIFS and GMOS spectra with
common TDE emission lines marked. With this deeper and higher resolution GMOS
spectrum we fit for the redshift using the methods of Hayden et al. (2021).
Using the weighted average of six binning choices from 5 to 20 Å we find a
redshift of $z=0.0734\pm 0.0004$ for the GMOS spectrum, which we will adopt
for the remainder of our analysis.

Figure 3: Optical spectra of ATLAS18mlw (black line) and the host galaxy WISEA
J073544.83+663717.3 (gray line). Vertical lines mark spectral features common
in TDEs and AGNs, with hydrogen lines in red, helium lines in blue, nitrogen
lines in green, and oxygen lines in orange. Both spectra have been binned in
$\sim$9 Å bins and the dichroic crossover region has been masked out.

##  3 Host-Galaxy Properties

Figure 4: Line measurements of the host galaxy WISEA J073544.83+663717.3
computed from the GMOS spectrum. Upper Left Panel: H$\alpha$ emission
equivalent width (EW), tracing ongoing star formation, as compared to the Lick
H$\delta_{A}$ absorption index, tracing star formation over the past Gyr. The
host WISEA J073544.83+663717.3 is represented by a red star and other TDE
hosts are shown with blue circles. WISEA J073544.83+663717.3 is a quiescent
Balmer strong galaxy, similar to many TDE hosts. Upper Right Panel: H$\alpha$
emission EW (WHα), as compared to log10([N II] / H$\alpha$), otherwise known
as the WHAN diagram (Cid Fernandes et al., 2011). Several lines separating
star-forming galaxies (SF), strong AGN (sAGN), weak AGN (wAGN), and passive
and “retired galaxies” (RG) are shown (Cid Fernandes et al., 2011). Lower Left
Panel: log10([O III] / H$\beta$) vs. log10([N II] / H$\alpha$) diagram
(Baldwin et al., 1981; Veilleux & Osterbrock, 1987). The solid line
theoretically separates AGNs (above) and H II-regions (below) (Kewley et al.,
2001). The dotted line empirically separates the same classifications
Kauffmann et al. (2003). Lower Right Panel: log10([O III] / H$\beta$) vs.
log10([S II] / H$\alpha$) diagram (Veilleux & Osterbrock, 1987). The solid
line represents the same separation as before (Kewley et al., 2001). The
diagonal dashed line separates Seyferts (above) and LINERs (below) (Kewley et
al., 2006). WISEA J073544.83+663717.3 has only upper limits on the H$\beta$
emission. Galaxies from the SDSS Data Release 8 (Eisenstein et al., 2011) are
shown in black to illustrate the broader distribution of galaxy properties in
these parameter spaces.

To better understand the host galaxy, we fit stellar population synthesis
models to the archival photometry of WISEA J073544.83+663717.3 (shown in Table
1), excluding the $FUV$ limit. We made use of the Fitting and Assessment of
Synthetic Templates (Fast; Kriek et al., 2009) code to obtain an SED for the
host galaxy. In this fit we assumed a Cardelli et al. (1988) extinction law
with $\text{R}_{\text{V}}=3.1$ and a foreground Galactic extinction of
$\text{A}_{\text{V}}=0.113$ mag (Schlafly & Finkbeiner, 2011), a Salpeter IMF
(Salpeter, 1955), an exponentially declining star-formation rate, and the
Bruzual & Charlot (2003) stellar population models. From our Fast fit, we find
that WISEA J073544.83+663717.3 has a stellar mass of
M${}_{*}=2.8^{+0.6}_{-0.4}$ $\times 10^{9}$ M⊙, an age of $2.8^{+0.7}_{-0.6}$
Gyr, and a star formation rate of SFR $=1.5^{+1.3}_{-1.1}\times 10^{-2}$ M⊙
yr-1.

From the scaling relationship between stellar mass and SMBH mass (McConnell &
Ma, 2013; Mendel et al., 2014), we estimate a central black hole mass of $\sim
10^{6.4}$ M⊙, where these scaling relations have a typical scatter of $\sim
0.4$ dex. This mass is close to the median BH mass on which TDEs occur ($\sim
10^{6}M_{\odot}$; Wevers et al., 2017; Ryu et al., 2020) and is also
consistent with the range of SMBH masses for which TDEs occur most efficiently
($\leq 10^{7}M_{\odot}$; van Velzen, 2018).

The GMOS spectrum provides additional insight into the host-galaxy properties.
Using the common scaling between SFR and H$\alpha$ luminosity from Osterbrock
& Ferland (2006), we obtain a SFR of $3.5\times 10^{-3}$ M⊙ yr-1, which is
smaller than, but consistent with, the value obtained from the SED fits.

In addition to the bulk host galaxy properties, we can constrain the presence
of AGN activity from the spectrum. We used the Penalized Pixel-Fitting code
(pPXF; Cappellari & Emsellem, 2004; Cappellari, 2017) to simultaneously fit a
stellar and gas component to the host galaxy spectrum, which allows for more
robust measurements of emission line fluxes. From pPXF we measure an [O III]
5007$\lambda$ emission line luminosity of $3.9\times 10^{39}$ erg s-1. This is
well below the median AGN [O III] line luminosity indicating that the host
galaxy does not harbor a strong AGN, but not ruling out the presence of low-
luminosity AGN (e.g., Zakamska et al., 2003; Bongiorno et al., 2010).

In addition to straightforward line-luminosity measurements, we make use of
well-accepted emission-line diagnostic diagrams to understand the host galaxy
properties and activity levels. Again using pPXF, we estimated emission line
fluxes and uncertainties from the GMOS spectrum. In the upper left panel of
Figure 4 we show the H$\alpha$ emission equivalent width (EW) as compared to
the Lick H$\delta_{A}$ absorption index. The H$\alpha$ EW traces current star
formation whereas the Lick H$\delta_{A}$ index traces star formation within
the past Gyr. This diagram is commonly used to select post-merger and post-
starburst systems based on their low current star formation rate but high
recent star formation (e.g., French et al., 2018). The host galaxy of
ATLAS18mlw has Lick H$\delta_{A}$ $\simeq 2.4\pm 1.0$ Å, as measured by
PyPHOT333https://mfouesneau.github.io/docs/pyphot/#, consistent with a
quiescent Balmer strong galaxy, similar to many other TDE host galaxies.

In the upper right panel of Figure 4 we show the H$\alpha$ emission EW (WHα),
as compared to log10([N II] / H$\alpha$), also known as the WHAN diagram (Cid
Fernandes et al., 2011). This diagram separates actively star forming and
active galaxies from passive galaxies using the H$\alpha$ emission EW. The
galaxies with high WHα are further separated by log10([N II] / H$\alpha$) into
pure star-forming galaxies and active galaxies. The host galaxy of ATLAS18mlw
resides in between the retired and passive galaxies regimes, suggesting that
the host does not harbor a strong AGN. Again, the host galaxy is similar in
terms of its location in the WHAN diagram to other TDE hosts.

In the bottom two panels of Figure 4 we show the traditional log10([N II] /
H$\alpha$) and log10([S II] / H$\alpha$) emission line ratio diagrams (Baldwin
et al., 1981; Veilleux & Osterbrock, 1987). H$\beta$ emission was not detected
in the host galaxy spectrum, and we therefore use 3$\sigma$ upper limits. In
both diagrams, the lower limit on log10([O III] / H$\beta$) lies in the
AGN/Seyfert region, but similar to other TDE hosts (Graur et al., 2018).
Addition, the large uncertainty on [S II] flux leads to ambiguity on the
classification of the source between a star-forming galaxy, a Seyfert, and a
LINER (low-ionization nuclear emission-line region; Heckman, 1980). Finally,
if we assume a typical Balmer decrement of 2.86 (Osterbrock, 1989) to estimate
the H$\beta$ flux from the observed H$\alpha$ flux, we obtain a line ratio of
log10([O III] / H$\beta$) $\simeq 0.9$.

Finally, we measured an X-ray flux at the position of the host galaxy using
archival ROSAT (ROentgen SATellite; Voges et al., 1999) images. We find a
3$\sigma$ upper limit of $<0.03$ counts s-1 in the $0.3-2.0$ keV bandpass.
Assuming an AGN-like photon index of 1.75 (Ricci et al., 2017) and the
Galactic column density along the line of sight of $3.63\times 10^{20}$ cm-2
(HI4PI Collaboration et al., 2016), we find a $0.3-10$ keV limit of
$<1.0\times 10^{-12}$ erg s-1 cm-2, which corresponds to a luminosity of
$<1.3\times 10^{43}$ erg s-1. This is fully consistent with a normal AGN, and
cannot solely rule out the presence of a strong AGN. Unfortunately, besides
ROSAT, there is little X-ray coverage of this host galaxy to further constrain
the prescence of an AGN. Nevertheless, combining emission line diagnostics and
the X-ray upper limit suggest that WISEA J073544.83+663717.3 is unlikely to
host a strong AGN.

##  4 Analysis

While there is relatively little data during the transient event, we can
nonetheless place constraints on the transient evolution and understand it in
the context of other TDEs. Fitting the H$\alpha$ emission as a single
Gaussian, we find a flux of $(7.3\pm 1.1)\times 10^{-15}$ erg cm-2 s-1
corresponding to a line luminosity of $(9.6\pm 1.4)\times 10^{40}$ erg s-1.
The equivalent width of the line relative to the nearby continuum is $112\pm
17$ Å, indicating a very strong emission line. The width of the spectral line
is $(1.5\pm 0.3)\times 10^{4}$ km s-1, all similar to known TDEs (e.g., van
Velzen et al., 2020; Gezari, 2021). Finally, because of the location of the
SNIFS dichroic crossover region at the redshift of ATLAS18mlw, it is not
possibly to cleanly detect and measure the properties of the tentative
H$\beta$ emission line.

A common analysis of optically-selected TDEs is blackbody fits to their
UV/optical SEDs to obtain effective temperatures and bolometric luminosities.
Unfortunately we do not have any UV photometry during the flare to probe the
high-energy emission from this event directly. It is known that TDEs without
UV photometry near their peaks can have sizable uncertainties on their SED
properties (e.g., Hinkle et al., 2020).

We can nonetheless make a rough estimate of these properties by fitting the
optical spectrum as a blackbody. To do this, we first corrected the SNIFS
spectrum for Galactic extinction. The GMOS host-galaxy spectrum does not cover
a wide enough range in wavelength to host-subtract the SNIFS transient
spectrum. Nevertheless, we fit the spectrum below $\sim 4900$ Å where the host
galaxy is faint and a large majority of the flux is from the transient.
Additionally, fitting below 4900 Å allowed us to avoid the affects of the
crossover region, where the SNIFS quick-look flux calibration is often poor
and artifacts can be present.

We fit the TDE spectrum with a Markov Chain Monte Carlo (MCMC) approach with
flat temperature priors of 1000 K $<$ T $<$ 55000 K. We assumed a blackbody
model for our fit and obtained the bolometric luminosity, radius, and
effective temperature. The $\chi^{2}$ per degree of freedom for this fit is
1.0, indicating a good fit to the blue portion of the spectrum. The spectrum
and corresponding blackbody fit are shown in Figure 5. The fitted luminosity
of L $=(1.7\pm 0.1)\times 10^{43}$ erg s-1, radius of R $=(1.1\pm 0.1)\times
10^{15}$ cm, and temperature of T $=(1.2\pm 0.1)\times 10^{4}$ K are also
shown on the figure. These uncertainties are very likely underestimated as
they only represent the statistical uncertainty on the fact and ignore the
systematic affects associated with fitting only optical data for hot sources.
As expected from the lack of host subtraction, the spectrum deviates from the
blackbody fit at redder wavelengths.

Figure 5: MCMC blackbody fit to the SNIFS spectrum of ATLAS18mlw. We present
the spectrum in $\lambda L_{\lambda}$ to show that the TDE is clearly
dominated by continuum emission in the blue and UV wavelengths. The best-
fitting luminosity, temperature, and radius are shown in the upper right
corner. The gray shaded area represents the portion of the spectrum that has
been fit as a blackbody. The gray dashed line represents 1 percent of the
Eddington luminosity for a 106.4 M⊙ SMBH. As before, the spectrum has been
binned in $\sim$9 Å bins and the dichroic crossover region has been masked
out.

Using the blackbody parameters computed above, we can better understand the
light curve of ATLAS18mlw in the context of other TDEs. Similar to previous
TDEs (Holoien et al., 2020; Hinkle et al., 2020, 2021a), we bolometrically
corrected the optical light curves shown in Fig. 2 using the fitted luminosity
at the time of the SNIFS spectrum. This allows us to compute an approximate
bolometric luminosity evolution for the TDE. This assumes a relatively flat
temperature evolution, which is often acceptable for a TDE (Hinkle et al.,
2021b; van Velzen et al., 2021; Gezari, 2021).

Following the procedure outlined in Hinkle et al. (2020), we computed the peak
luminosity and decline rate measured at 40 days for ATLAS18mlw. Given the
small number of epochs near peak light, we used a Monte Carlo approach to
estimate the peak luminosity, fitting the data near peak 5000 times, each time
with the luminosities perturbed by their uncertainties assuming Gaussian
errors. This yielded a peak luminosity of log(L [erg s-1]) = $43.5\pm 0.2$. To
ensure a robust uncertainty estimate given the lack of UV data, this error
represents the 90% confidence interval of the luminosities computed from our
Monte Carlo estimation. Estimating a time of peak light with the same
procedure yields MJD${}_{peak}=58204.4\pm 6.8$, consistent with the
approximate optical peak seen in Fig. 2. Given the SMBH mass estimated from
the host galaxy SED fitting, this represents an Eddington ratio of $\sim 10$%,
similar to other TDEs (e.g., Wevers et al., 2017). We next computed the
decline rate and found $\Delta L_{40}=-0.7\pm 0.2$. The peak luminosity and
decline rate of ATLAS18mlw, along with well-studied TDEs and ambiguous nuclear
transients (ANTs) are shown in Figure 6, with ATLAS18mlw broadly following the
trend seen for TDEs (e.g., Hinkle et al., 2020, 2021a).

Figure 6: Peak bolometric luminosity versus the decline rate of several TDEs
(blue squares) and ANTs (gray circles). The decline rate $\Delta$L40 is the
difference between the log of the peak luminosity and the log of the
luminosity at 40 days after peak (Hinkle et al., 2020). ATLAS18mlw is shown
with a red star, and the uncertainties here are 90% confidence intervals.

##  5 Discussion

The observational data of ATLAS18mlw appears similar to many other TDEs, with
no significant archival host variability, a likely post-merger host galaxy, a
strong blue continuum with broad emission lines near peak, and a timescale of
a few months in the optical. While we do not have direct measurements of the
UV flux, the presence of strong H$\alpha$ emission indicates that there must
be significant far UV flux in the neighborhood of the transient, typical of
TDEs (e.g., Roth & Kasen, 2018) and Seyfert galaxies (e.g., Antonucci, 1993).
In fact, the line width and luminosity of the H$\alpha$ line seen for
ATLAS18mlw is similar to other TDEs near peak like iPTF16fnl (Blagorodnova et
al., 2017), iPTF16axa (Hung et al., 2017), PS18kh (Holoien et al., 2019), and
ASASSN-19dj (Hinkle et al., 2021b). Additionally, the appearance of the broad
H$\alpha$ line at an approximate phase of nine days before peak is among the
earliest detected. Many TDEs show an increase in line luminosity near or
shortly after peak light before declining (Charalampopoulos et al., 2021),
although here we only have a single spectrum with which to constrain the line
emission.

Transients associated with AGN accretion rate changes also commonly exhibit
broad H$\alpha$ emission (e.g., Gezari et al., 2017; Neustadt et al., 2020;
Frederick et al., 2020). As such, we must consider the possibility that
ATLAS18mlw is a smooth AGN flare. Firstly, the short overall timescale of
ATLAS18mlw is less than the roughly year-long timescales of AGN-related flares
(Hinkle et al., 2021a). Additionally, the spectrum of ATLAS18mlw is well-fit
as blackbody, consistent with TDEs (e.g., Holoien et al., 2014a, 2016),
whereas most AGN are better fit by power-laws (e.g., Vanden Berk et al., 2001;
Neustadt et al., 2020). Finally, the lack of significant optical variability
and weak H$\alpha$ emission from the host galaxy are inconsistent with the
presence of an AGN.

Using the optical spectroscopic classification scheme of Leloudas et al.
(2019) and van Velzen et al. (2021), this TDE would likely be a TDE-H due to
the dominance of H$\alpha$ over any other spectral feature. TDE-H objects tend
to have broader lines than TDE-Bowen objects (Charalampopoulos et al., 2021),
consistent with the H$\alpha$ line width computed for ATLAS18mlw. There does
appear to be a He II line, but it is narrow and seen in the host galaxy
spectrum. However, some TDEs, like ASASSN-15oi (Holoien et al., 2018) and
ASASSN-19dj (Hinkle et al., 2021b) show signs for changes in spectral features
and the properties of observed emission lines, which we cannot probe for
ATLAS18mlw with our single spectrum.

The temperature derived from our blackbody fit, log(T[K]) = $4.06\pm 0.01$, is
towards the low end of temperatures for the TDE-H subclass (van Velzen et al.,
2021). However, the lack of UV photometry makes the systematic errors on our
temperature significantly higher than the statistical errors from the
blackbody fit alone. Nonetheless, both the TDE-He and TDE-Bowen classes have
much higher median temperatures, in support of this object belonging to the
TDE-H class. Additionally, the estimated blackbody radius of log(R [cm]) =
$15.07\pm 0.02$ is fully consistent with the radii of the TDE-H objects. van
Velzen et al. (2021) find that the separation between TDE-H and TDE-Bowen is
strongest for the blackbody radius, again supporting a TDE-H nature. Still, if
our crude blackbody fit without UV photometry underestimates the temperature,
then the blackbody radius may be overestimated for a similar luminosity.

In terms of the bolometric evolution, we find that ATLAS18mlw is most similar
to the “faint and fast” TDEs iPTF16fnl (Blagorodnova et al., 2017; Brown et
al., 2018) and ZTF19abzrhgq (AT2019qiz; Nicholl et al., 2020). Of the three,
ATLAS18mlw is the lowest luminosity, although the bolometric luminosity
estimate is the weakest as there is no UV data to confirm the temperature.
ATLAS18mlw is also the only likely TDE-H object of these three TDEs. The
decline rate of ATLAS18mlw indicates that it is a fast evolving TDE, with only
three sources declining faster. These are iPTF16fnl, ZTF19abzrhgq, and
ASASSN-14ae (Holoien et al., 2014b). It is also interesting to note that the
position of ATLAS18mlw in peak-luminosity/decline-rate space is much more
similar to the TDEs than any of the ambiguous nuclear transients (ANTs) shown
on Fig. 6. This result supports using the log(Lpeak) – $\Delta$L40 parameter
space to help us differentiate between different classes of transients even
without abundant data during the event, as will be the case for many
discoveries by the Legacy Survey of Space and Time (LSST; Ivezić et al., 2019)
on the Vera Rubin Observatory.

For the case of ATLAS18mlw, the late identification as a TDE precluded any
follow-up observations, highlighting the need for rapid response to newly-
discovered nuclear transients. As the discovery of TDEs continues to expand,
it is important to prioritize data that best helps us understand physics
driving the observed emission. Optical spectroscopy and high-energy UV/X-ray
photometry, ideally beginning prior to peak light, are necessary to understand
the detailed properties of individual events.

This goal requires the expansion of spectroscopic classification surveys to
match the growing flood of transient discoveries each night. Even now, the
number of discovered transients is roughly an order of magnitude higher than
the number of transients classified spectroscopically. One obvious improvement
is the robotizing of existing classification surveys to improve efficiency.
For instance, the SCAT survey is currently classically scheduled and observed
with 323 classifications in the last year, representing $\sim 13$% of all
classified transients. After the ongoing robotizing of the UH88
telescope444This effort is funded by an NSF Major Research Instrumentation
Program (MRI) grant (Award Number: 1920392). is completed in late 2022, this
number is expected to increase to $\sim 1200$ per year.

Sources like ATLAS18mlw show that even with minimal data, important
conclusions can be drawn and comparisons to very well-studied objects can be
made. In particular, the host galaxies of sources like ATLAS18mlw can and
should be incorporated into TDE population studies as these samples grow.

We thank Alexa Anderson for helpful comments on the manuscript. J.T.H. and
this work was supported by NASA award 80NSSC21K0136. M.A.T. acknowledges
support from the DOE CSGF through grant DE-SC0019323. Support for T.W.-S.H.
was provided by NASA through the NASA Hubble Fellowship grant HST-
HF2-51458.001-A awarded by the Space Telescope Science Institute, which is
operated by the Association of Universities for Research in Astronomy, Inc.,
for NASA, under contract NAS5-26555. B.J.S is supported by NSF grants
AST-1908952, AST-1920392, AST-1911074, and NASA award 80NSSC19K1717. Support
for G.A. was provided by the Director, Office of Science, Office of High
Energy Physics of the U.S. Department of Energy under Contract No. DE-
AC025CH11231 This work uses data from the University of Hawaii’s ATLAS
project, funded through NASA grants NN12AR55G, 80NSSC18K0284, and
80NSSC18K1575, with contributions from the Queen’s University Belfast, STScI,
the South African Astronomical Observatory, and the Millennium Institute of
Astrophysics, Chile. This paper includes data collected by the TESS mission,
which are publicly available from the Mikulski Archive for Space Telescopes
(MAST). Funding for the TESS mission is provided by NASA’s Science Mission
directorate. Parts of this research were supported by the Australian Research
Council Centre of Excellence for All Sky Astrophysics in 3 Dimensions (ASTRO
3D), through project number CE170100013. This work is based on observations
made by ATLAS and Gemini. The authors wish to recognize and acknowledge the
very significant cultural role and reverence that the summits of Haleakalā and
Maunakea have always had within the indigenous Hawaiian community. We are most
fortunate to have the opportunity to conduct observations from these
mountains.

## References

  * Aguado et al. (2019) Aguado, D. S., Ahumada, R., Almeida, A., et al. 2019, ApJS, 240, 23, doi: 10.3847/1538-4365/aaf651
  * Alard (2000) Alard, C. 2000, AAPS, 144, 363, doi: 10.1051/aas:2000214
  * Alard & Lupton (1998) Alard, C., & Lupton, R. H. 1998, ApJ, 503, 325, doi: 10.1086/305984
  * Aldering et al. (2006) Aldering, G., Antilogus, P., Bailey, S., et al. 2006, ApJ, 650, 510, doi: 10.1086/507020
  * Antonucci (1993) Antonucci, R. 1993, ARA&A, 31, 473, doi: 10.1146/annurev.aa.31.090193.002353
  * Arcavi et al. (2014) Arcavi, I., Gal-Yam, A., Sullivan, M., et al. 2014, ApJ, 793, 38, doi: 10.1088/0004-637X/793/1/38
  * Auchettl et al. (2017) Auchettl, K., Guillochon, J., & Ramirez-Ruiz, E. 2017, ApJ, 838, 149, doi: 10.3847/1538-4357/aa633b
  * Auchettl et al. (2018) Auchettl, K., Ramirez-Ruiz, E., & Guillochon, J. 2018, ApJ, 852, 37, doi: 10.3847/1538-4357/aa9b7c
  * Baldwin et al. (1981) Baldwin, J. A., Phillips, M. M., & Terlevich, R. 1981, PASP, 93, 5, doi: 10.1086/130766
  * Barbarino et al. (2019) Barbarino, C., Carracedo, A. S., Tartaglia, L., & Yaron, O. 2019, Transient Name Server Classification Report, 2019-287, 1
  * Barron et al. (2008) Barron, J. T., Stumm, C., Hogg, D. W., Lang, D., & Roweis, S. 2008, AJ, 135, 414, doi: 10.1088/0004-6256/135/1/414
  * Bellm et al. (2019) Bellm, E. C., Kulkarni, S. R., Graham, M. J., et al. 2019, PASP, 131, 018002, doi: 10.1088/1538-3873/aaecbe
  * Bennett et al. (2014) Bennett, C. L., Larson, D., Weiland, J. L., & Hinshaw, G. 2014, ApJ, 794, 135, doi: 10.1088/0004-637X/794/2/135
  * Blagorodnova et al. (2017) Blagorodnova, N., Gezari, S., Hung, T., et al. 2017, ApJ, 844, 46, doi: 10.3847/1538-4357/aa7579
  * Bongiorno et al. (2010) Bongiorno, A., Mignoli, M., Zamorani, G., et al. 2010, A&A, 510, A56, doi: 10.1051/0004-6361/200913229
  * Brown et al. (2017) Brown, J. S., Holoien, T. W.-S., Auchettl, K., et al. 2017, MNRAS, 466, 4904, doi: 10.1093/mnras/stx033
  * Brown et al. (2018) Brown, J. S., Kochanek, C. S., Holoien, T. W.-S., et al. 2018, MNRAS, 473, 1130, doi: 10.1093/mnras/stx2372
  * Bruzual & Charlot (2003) Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000, doi: 10.1046/j.1365-8711.2003.06897.x
  * Cappellari (2017) Cappellari, M. 2017, MNRAS, 466, 798, doi: 10.1093/mnras/stw3020
  * Cappellari & Emsellem (2004) Cappellari, M., & Emsellem, E. 2004, PASP, 116, 138, doi: 10.1086/381875
  * Cardelli et al. (1988) Cardelli, J. A., Clayton, G. C., & Mathis, J. S. 1988, ApJ, 329, L33, doi: 10.1086/185171
  * Chambers et al. (2016) Chambers, K. C., Magnier, E. A., Metcalfe, N., et al. 2016, ArXiv e-prints.  https://arxiv.org/abs/1612.05560
  * Charalampopoulos et al. (2021) Charalampopoulos, P., Leloudas, G., Malesani, D. B., et al. 2021, arXiv e-prints, arXiv:2109.00016.  https://arxiv.org/abs/2109.00016
  * Cid Fernandes et al. (2011) Cid Fernandes, R., Stasińska, G., Mateus, A., & Vale Asari, N. 2011, MNRAS, 413, 1687, doi: 10.1111/j.1365-2966.2011.18244.x
  * Eisenstein et al. (2011) Eisenstein, D. J., Weinberg, D. H., Agol, E., et al. 2011, AJ, 142, 72, doi: 10.1088/0004-6256/142/3/72
  * Evans & Kochanek (1989) Evans, C. R., & Kochanek, C. S. 1989, ApJ, 346, L13, doi: 10.1086/185567
  * Frederick et al. (2020) Frederick, S., Gezari, S., Graham, M. J., et al. 2020, arXiv e-prints, arXiv:2010.08554.  https://arxiv.org/abs/2010.08554
  * Fremling et al. (2020) Fremling, C., Miller, A. A., Sharma, Y., et al. 2020, ApJ, 895, 32, doi: 10.3847/1538-4357/ab8943
  * Fremling et al. (2021) Fremling, C., Hall, X. J., Coughlin, M. W., et al. 2021, ApJ, 917, L2, doi: 10.3847/2041-8213/ac116f
  * French (2021) French, K. D. 2021, PASP, 133, 072001, doi: 10.1088/1538-3873/ac0a59
  * French et al. (2016) French, K. D., Arcavi, I., & Zabludoff, A. 2016, ApJ, 818, L21, doi: 10.3847/2041-8205/818/1/L21
  * French et al. (2020a) French, K. D., Arcavi, I., Zabludoff, A. I., et al. 2020a, ApJ, 891, 93, doi: 10.3847/1538-4357/ab7450
  * French et al. (2020b) French, K. D., Wevers, T., Law-Smith, J., Graur, O., & Zabludoff, A. I. 2020b, Space Sci. Rev., 216, 32, doi: 10.1007/s11214-020-00657-y
  * French et al. (2018) French, K. D., Yang, Y., Zabludoff, A. I., & Tremonti, C. A. 2018, ApJ, 862, 2, doi: 10.3847/1538-4357/aacb2d
  * Gafton & Rosswog (2019) Gafton, E., & Rosswog, S. 2019, MNRAS, 487, 4790, doi: 10.1093/mnras/stz1530
  * Gezari (2021) Gezari, S. 2021, arXiv e-prints, arXiv:2104.14580.  https://arxiv.org/abs/2104.14580
  * Gezari et al. (2017) Gezari, S., Hung, T., Cenko, S. B., et al. 2017, ApJ, 835, 144, doi: 10.3847/1538-4357/835/2/144
  * Graham et al. (2001) Graham, A. W., Erwin, P., Caon, N., & Trujillo, I. 2001, ApJ, 563, L11, doi: 10.1086/338500
  * Graur et al. (2018) Graur, O., French, K. D., Zahid, H. J., et al. 2018, ApJ, 853, 39, doi: 10.3847/1538-4357/aaa3fd
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357, doi: 10.1038/s41586-020-2649-2
  * Hayden et al. (2021) Hayden, B., Rubin, D., Boone, K., et al. 2021, ApJ, 912, 87, doi: 10.3847/1538-4357/abed4d
  * Heckman (1980) Heckman, T. M. 1980, A&A, 500, 187
  * HI4PI Collaboration et al. (2016) HI4PI Collaboration, Ben Bekhti, N., Flöer, L., et al. 2016, A&A, 594, A116, doi: 10.1051/0004-6361/201629178
  * Hinkle et al. (2021a) Hinkle, J. T., Holoien, T. W. S., Shappee, B. J., & Auchettl, K. 2021a, ApJ, 910, 83, doi: 10.3847/1538-4357/abe4d8
  * Hinkle et al. (2020) Hinkle, J. T., Holoien, T. W. S., Shappee, B. J., et al. 2020, ApJ, 894, L10, doi: 10.3847/2041-8213/ab89a2
  * Hinkle et al. (2021b) Hinkle, J. T., Holoien, T. W. S., Auchettl, K., et al. 2021b, MNRAS, 500, 1673, doi: 10.1093/mnras/staa3170
  * Ho (2008) Ho, L. C. 2008, ARA&A, 46, 475, doi: 10.1146/annurev.astro.45.051806.110546
  * Holoien et al. (2018) Holoien, T. W.-S., Brown, J. S., Auchettl, K., et al. 2018, MNRAS, 480, 5689, doi: 10.1093/mnras/sty2273
  * Holoien et al. (2014a) Holoien, T. W.-S., Prieto, J. L., Stanek, K. Z., et al. 2014a, ApJ, 785, L35, doi: 10.1088/2041-8205/785/2/L35
  * Holoien et al. (2014b) Holoien, T. W.-S., Prieto, J. L., Bersier, D., et al. 2014b, MNRAS, 445, 3263, doi: 10.1093/mnras/stu1922
  * Holoien et al. (2016) Holoien, T. W.-S., Kochanek, C. S., Prieto, J. L., et al. 2016, MNRAS, 455, 2918, doi: 10.1093/mnras/stv2486
  * Holoien et al. (2019) Holoien, T. W. S., Huber, M. E., Shappee, B. J., et al. 2019, ApJ, 880, 120, doi: 10.3847/1538-4357/ab2ae1
  * Holoien et al. (2020) Holoien, T. W. S., Auchettl, K., Tucker, M. A., et al. 2020, arXiv e-prints, arXiv:2003.13693.  https://arxiv.org/abs/2003.13693
  * Hook et al. (2004) Hook, I. M., Jørgensen, I., Allington-Smith, J. R., et al. 2004, PASP, 116, 425, doi: 10.1086/383624
  * Hung et al. (2017) Hung, T., Gezari, S., Blagorodnova, N., et al. 2017, ApJ, 842, 29, doi: 10.3847/1538-4357/aa7337
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science & Engineering, 9, 90, doi: 10.1109/MCSE.2007.55
  * Ivezić et al. (2019) Ivezić, Ž., Kahn, S. M., Tyson, J. A., et al. 2019, ApJ, 873, 111, doi: 10.3847/1538-4357/ab042c
  * Kauffmann et al. (2003) Kauffmann, G., Heckman, T. M., Tremonti, C., et al. 2003, MNRAS, 346, 1055, doi: 10.1111/j.1365-2966.2003.07154.x
  * Kelly (2007) Kelly, B. C. 2007, ApJ, 665, 1489, doi: 10.1086/519947
  * Kewley et al. (2001) Kewley, L. J., Dopita, M. A., Sutherland, R. S., Heisler, C. A., & Trevena, J. 2001, ApJ, 556, 121, doi: 10.1086/321545
  * Kewley et al. (2006) Kewley, L. J., Groves, B., Kauffmann, G., & Heckman, T. 2006, MNRAS, 372, 961, doi: 10.1111/j.1365-2966.2006.10859.x
  * Kochanek et al. (2017) Kochanek, C. S., Shappee, B. J., Stanek, K. Z., et al. 2017, PASP, 129, 104502, doi: 10.1088/1538-3873/aa80d9
  * Kriek et al. (2009) Kriek, M., van Dokkum, P. G., Labbé, I., et al. 2009, ApJ, 700, 221, doi: 10.1088/0004-637X/700/1/221
  * Lacerda et al. (2020) Lacerda, E. A. D., Sánchez, S. F., Cid Fernandes, R., et al. 2020, MNRAS, 492, 3073, doi: 10.1093/mnras/staa008
  * Lang et al. (2010) Lang, D., Hogg, D. W., Mierle, K., Blanton, M., & Roweis, S. 2010, AJ, 139, 1782, doi: 10.1088/0004-6256/139/5/1782
  * Lantz et al. (2004) Lantz, B., Aldering, G., Antilogus, P., et al. 2004, in Proc. SPIE, Vol. 5249, Optical Design and Engineering, ed. L. Mazuray, P. J. Rogers, & R. Wartmann, 146–155, doi: 10.1117/12.512493
  * Law-Smith et al. (2017) Law-Smith, J., Ramirez-Ruiz, E., Ellison, S. L., & Foley, R. J. 2017, ApJ, 850, 22, doi: 10.3847/1538-4357/aa94c7
  * Leloudas et al. (2019) Leloudas, G., Dai, L., Arcavi, I., et al. 2019, arXiv e-prints.  https://arxiv.org/abs/1903.03120
  * Martin et al. (2005) Martin, D. C., Fanson, J., Schiminovich, D., et al. 2005, ApJ, 619, L1, doi: 10.1086/426387
  * McConnell & Ma (2013) McConnell, N. J., & Ma, C.-P. 2013, ApJ, 764, 184, doi: 10.1088/0004-637X/764/2/184
  * Mendel et al. (2014) Mendel, J. T., Simard, L., Palmer, M., Ellison, S. L., & Patton, D. R. 2014, ApJS, 210, 3, doi: 10.1088/0067-0049/210/1/3
  * Million et al. (2016) Million, C., Fleming, S. W., Shiao, B., et al. 2016, ApJ, 833, 292, doi: 10.3847/1538-4357/833/2/292
  * Mockler et al. (2019) Mockler, B., Guillochon, J., & Ramirez-Ruiz, E. 2019, ApJ, 872, 151, doi: 10.3847/1538-4357/ab010f
  * Mockler et al. (2022) Mockler, B., Twum, A. A., Auchettl, K., et al. 2022, ApJ, 924, 70, doi: 10.3847/1538-4357/ac35d5
  * Neustadt et al. (2020) Neustadt, J. M. M., Holoien, T. W. S., Kochanek, C. S., et al. 2020, MNRAS, 494, 2538, doi: 10.1093/mnras/staa859
  * Nicholl et al. (2020) Nicholl, M., Wevers, T., Oates, S. R., et al. 2020, arXiv e-prints, arXiv:2006.02454.  https://arxiv.org/abs/2006.02454
  * Osterbrock (1989) Osterbrock, D. E. 1989, Astrophysics of Gaseous Nebulae and Active Galactic Nuclei
  * Osterbrock & Ferland (2006) Osterbrock, D. E., & Ferland, G. J. 2006, Astrophysics of gaseous nebulae and active galactic nuclei
  * Phinney (1989) Phinney, E. S. 1989, Nature, 340, 595, doi: 10.1038/340595a0
  * Prieto et al. (2016) Prieto, J. L., Krühler, T., Anderson, J. P., et al. 2016, ApJ, 830, L32, doi: 10.3847/2041-8205/830/2/L32
  * Rees (1988) Rees, M. J. 1988, Nature, 333, 523, doi: 10.1038/333523a0
  * Ricci et al. (2017) Ricci, C., Trakhtenbrot, B., Koss, M. J., et al. 2017, ApJS, 233, 17, doi: 10.3847/1538-4365/aa96ad
  * Ricker et al. (2015) Ricker, G. R., Winn, J. N., Vanderspek, R., et al. 2015, Journal of Astronomical Telescopes, Instruments, and Systems, 1, 014003, doi: 10.1117/1.JATIS.1.1.014003
  * Roth & Kasen (2018) Roth, N., & Kasen, D. 2018, ApJ, 855, 54, doi: 10.3847/1538-4357/aaaec6
  * Ryu et al. (2020) Ryu, T., Krolik, J., & Piran, T. 2020, arXiv e-prints, arXiv:2007.13765.  https://arxiv.org/abs/2007.13765
  * Salpeter (1955) Salpeter, E. E. 1955, ApJ, 121, 161, doi: 10.1086/145971
  * Schlafly & Finkbeiner (2011) Schlafly, E. F., & Finkbeiner, D. P. 2011, ApJ, 737, 103, doi: 10.1088/0004-637X/737/2/103
  * Shappee et al. (2014) Shappee, B. J., Prieto, J. L., Grupe, D., et al. 2014, ApJ, 788, 48, doi: 10.1088/0004-637X/788/1/48
  * Smartt et al. (2015) Smartt, S. J., Valenti, S., Fraser, M., et al. 2015, A&A, 579, A40, doi: 10.1051/0004-6361/201425237
  * Smith et al. (2020) Smith, K. W., Smartt, S. J., Young, D. R., et al. 2020, arXiv e-prints, arXiv:2003.09052.  https://arxiv.org/abs/2003.09052
  * Stone & Metzger (2016) Stone, N. C., & Metzger, B. D. 2016, MNRAS, 455, 859, doi: 10.1093/mnras/stv2281
  * Stone & van Velzen (2016) Stone, N. C., & van Velzen, S. 2016, ApJ, 825, L14, doi: 10.3847/2041-8205/825/1/L14
  * Tonry et al. (2018) Tonry, J. L., Denneau, L., Heinze, A. N., et al. 2018, PASP, 130, 064505, doi: 10.1088/1538-3873/aabadf
  * Tucker et al. (2018) Tucker, M. A., Huber, M., Shappee, B. J., et al. 2018, The Astronomer’s Telegram, 11444, 1
  * Ulmer (1999) Ulmer, A. 1999, ApJ, 514, 180, doi: 10.1086/306909
  * Vallely et al. (2021) Vallely, P. J., Kochanek, C. S., Stanek, K. Z., Fausnaugh, M., & Shappee, B. J. 2021, MNRAS, 500, 5639, doi: 10.1093/mnras/staa3675
  * Vallely et al. (2019) Vallely, P. J., Fausnaugh, M., Jha, S. W., et al. 2019, arXiv e-prints.  https://arxiv.org/abs/1903.08665
  * van Velzen (2018) van Velzen, S. 2018, ApJ, 852, 72, doi: 10.3847/1538-4357/aa998e
  * van Velzen & Farrar (2014) van Velzen, S., & Farrar, G. R. 2014, ApJ, 792, 53, doi: 10.1088/0004-637X/792/1/53
  * van Velzen et al. (2020) van Velzen, S., Holoien, T. W. S., Onori, F., Hung, T., & Arcavi, I. 2020, Space Sci. Rev., 216, 124, doi: 10.1007/s11214-020-00753-z
  * van Velzen et al. (2019) van Velzen, S., Stone, N. C., Metzger, B. D., et al. 2019, ApJ, 878, 82, doi: 10.3847/1538-4357/ab1844
  * van Velzen et al. (2021) van Velzen, S., Gezari, S., Hammerstein, E., et al. 2021, ApJ, 908, 4, doi: 10.3847/1538-4357/abc258
  * Vanden Berk et al. (2001) Vanden Berk, D. E., Richards, G. T., Bauer, A., et al. 2001, AJ, 122, 549, doi: 10.1086/321167
  * Vanderspek et al. (2018) Vanderspek, R., Doty, J., Fausnaugh, M., et al. 2018, TESS Instrument Handbook, Tech. rep., Kavli Institute for Astrophysics and Space Science, Massachusetts Institute of Technology
  * Veilleux & Osterbrock (1987) Veilleux, S., & Osterbrock, D. E. 1987, ApJS, 63, 295, doi: 10.1086/191166
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261, doi: 10.1038/s41592-019-0686-2
  * Voges et al. (1999) Voges, W., Aschenbach, B., Boller, T., et al. 1999, A&A, 349, 389
  * Wang & Merritt (2004) Wang, J., & Merritt, D. 2004, ApJ, 600, 149, doi: 10.1086/379767
  * Wevers et al. (2017) Wevers, T., van Velzen, S., Jonker, P. G., et al. 2017, MNRAS, 471, 1694, doi: 10.1093/mnras/stx1703
  * Wright (2006) Wright, E. L. 2006, PASP, 118, 1711, doi: 10.1086/510102
  * Wright et al. (2010) Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868, doi: 10.1088/0004-6256/140/6/1868
  * Wyrzykowski et al. (2012) Wyrzykowski, Ł., Hodgkin, S., Blogorodnova, N., Koposov, S., & Burgon, R. 2012, in 2nd Gaia Follow-up Network for Solar System Objects, 21\.  https://arxiv.org/abs/1210.5007
  * Zakamska et al. (2003) Zakamska, N. L., Strauss, M. A., Krolik, J. H., et al. 2003, AJ, 126, 2125, doi: 10.1086/378610


# Predicting transitions in cooperation levels from network connectivity

A. Zhuk1,+, I. Sendiña-Nadal2,3,+, I. Leyva2,3,+, D. Musatov1,4,5, A.M.
Raigorodskii1,5,6,7, M. Perc8,9,10,11, and S. Boccaletti1,2,12 1Moscow
Institute of Physics and Technology (National Research University), 9
Institutskiy per., Dolgoprudny, Moscow Region, 141701, Russia  2Universidad
Rey Juan Carlos, Calle Tulipán s/n, 28933 Móstoles, Madrid, Spain  3Center for
Biomedical Technology, Universidad Politécnica de Madrid, 28223 Pozuelo de
Alarcón, Madrid, Spain  4Russian Academy of National Economy and Public
Administration, pr. Vernadskogo, 84, Moscow, 119606, Russia  5 Caucasus
Mathematical Center at Adyghe State University, ul. Pervomaiskaya, 208,
Maykop, 385000, Russia  6Mechanics and Mathematics Faculty, Moscow State
University, Leninskie Gory, 1, Moscow, 119991, Russia  7Institute of
Mathematics and Computer Science, Buryat State University, ul. Ranzhurova, 5,
Ulan-Ude, 670000, Russia  8Faculty of Natural Sciences and Mathematics,
University of Maribor, Koroška cesta 160, 2000 Maribor, Slovenia  9Department
of Medical Research, China Medical University Hospital, China Medical
University, 404332 Taichung, Taiwan  10Alma Mater Europaea, Slovenska ulica
17, 2000 Maribor, Slovenia  11Complexity Science Hub Vienna,
Josefstädterstraße 39, 1080 Vienna, Austria  12CNR - Institute of Complex
Systems, Via Madonna del Piano 10, I-50019 Sesto Fiorentino, Italy
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Networks determine our social circles and the way we cooperate with others. We
know that topological features like hubs and degree assortativity affect
cooperation, and we know that cooperation is favoured if the benefit of the
altruistic act divided by the cost exceeds the average number of neighbours.
However, a simple rule that would predict cooperation transitions on an
arbitrary network has not yet been presented. Here we show that the unique
sequence of degrees in a network can be used to predict at which game
parameters major shifts in the level of cooperation can be expected, including
phase transitions from absorbing to mixed strategy phases. We use the
evolutionary prisoner’s dilemma game on random and scale-free networks to
demonstrate the prediction, as well as its limitations and possible pitfalls.
We observe good agreements between the predictions and the results obtained
with concurrent and Monte Carlo methods for the update of the strategies, thus
providing a simple and fast way to estimate the outcome of evolutionary social
dilemmas on arbitrary networks without the need of actually playing the game.

††: New J. Phys.

,
+These Authors contributed equally

##  1 Introduction

In 1992 Nowak and May observed that cooperators form compact clusters and can
thereby withstand invading defectors in an iterated prisoner’s dilemma game on
a square lattice [1]. This was a fascinating discovery because cooperators
should have died out in agreement with the Nash equilibrium of the game [2,
3]. This would, in fact, always happen in a well-mixed population, but not on
a network [4, 5, 6]. The fact that the structure of a network can positively
affect the evolution of cooperation is today known as network reciprocity [1,
7], and it constitutes one of five key mechanisms for cooperators to survive
social dilemmas [8].

With the advent of network science at the turn of the 21st century [9, 10,
11], the field of evolutionary games on networks rapidly gained on popularity,
and various complex networks have been studied for their impact on the
evolution of cooperation, including scale-free [12, 13, 14, 15, 16, 17, 18,
19, 20, 21, 22, 23], small-world [24, 25, 26, 27, 28, 29, 30], hierarchical
[31, 32], coevolving [33, 34, 35, 36, 37, 38], and empirical social networks
[39, 40, 41]. Later on, interdependent and multilayer networks have emerged as
an important new paradigm in network science [42, 43], and these were also
considered prominently in cooperation research [44, 45, 46, 47, 48, 49, 50,
51, 52], including the discovery of interdependent network reciprocity as an
extension of the traditional case [53]. Most recently, the focus has been
shifting once more, this time to higher-order networks [54], where
evolutionary dynamics has been considered [55, 56, 57].

However, despite the wealth of research concerning the evolution of
cooperation on networks, fundamental results have been relatively scarce.
Tanimoto [58] and Wang et al. [59] introduced universal scaling parameters to
account for the strength of a dilemma disrupting the promotion of cooperation
that can evaluate the outcome of a game independently of the reciprocity
mechanisms, and showed how these mechanisms relax the dilemma strength [60].
Ohtsuki et al. [61], for instance, presented a simple rule for the evolution
of cooperation on graphs and social networks, proving that natural selection
favours cooperation if the benefit of the altruistic act, divided by the cost,
exceeds the average number of neighbours. Moreover, Allen et al. [62] derived
a general formula for pairwise games under weak selection that applies to any
graph or social network.

In this paper, we derive a simple conjecture for cooperation transitions in
evolutionary social dilemmas based on the unique sequence of degrees in a
network. As we will show, the conjecture predicts well at which game
parameters important shifts in the level of cooperation can be expected,
including phase transitions from absorbing to mixed strategy phases [63, 64].

In what follows, we first present the derivation of the conjecture and then
proceed with showing the results obtained for the evolutionary prisoner’s
dilemma game on random and scale-free networks. We conclude by discussing the
implications of our research to optimise large-scale simulations of
evolutionary processes on complex networks and the possibilities for the
generalisation of the conjecture to related subjects.

##  2 Game and conjecture formulation

We start by considering a population of $N$ agents playing the Prisoner’s
Dilemma (PD) game [65]. At each step of the game, agents are allowed to play
different strategies, that is, a node can cooperate with some neighbours and
defect with others [66], also known as the mixed strategy system with
different offers [67, 68].

Connections between players follow the network structure prescribed by the
adjacency matrix $A=(a_{uv})$, where $a_{uv}=1$ if players $u$ and $v$ are
linked, and $a_{uv}=0$ otherwise.

The payoffs in the PD game are extracted from the following matrix:

${PD}(b,c)=\begin{array}[]{c|cc}&{\rm C}&{\rm D}\\\ \hline\cr{\rm C}&b-c&-c\\\
\hline\cr{\rm D}&b&0\\\ \end{array}$ (1)

which means that a cooperator (C) altruistically loses an amount $c>0$ when
playing with a defector (D), who wins an amount $b>c$ in such interaction. An
interplay between two cooperators results in a mutual benefit of $b-c$,
whereas playing between two defectors results in a null interchange. Note that
a defector never looses.

Without loss of generality, one can take $c=1$ ($b>1$) and re-scale the payoff
matrix as:

$P(b)={PD}(b,1)/b=\begin{pmatrix}1-1/b&-1/b\\\ 1&0\\\ \end{pmatrix}.$ (2)

According to the universal scaling for the dilemma strength introduced in Ref.
[59], $D^{\prime}_{g}=D^{\prime}_{r}=1/(b-1)$, the parameter $b$ measures the
dilemma weakness.

For interacting players $u$ and $v$, one can define the strategy variable
$x_{uv}$ as

$x_{uv}=\left\\{\begin{array}[]{lll}1&\mbox{if}&u\;\mbox{cooperates with }v\\\
0&\mbox{if}&u\;\mbox{defects }v\end{array}\right.$

and the strategy vector as $\overline{x}_{uv}=\big{(}x_{uv},1-x_{uv}\big{)}$.
In each round of the game, the payoff $g_{uv}$ earned by each player $u$ from
its interaction with $v$ is:

$g_{uv}=\overline{x}_{uv}P(b)\overline{x}_{vu}^{T}.$ (3)

The average payoff $g_{u}$ collected by player $u$ after an interaction round
with its $k_{u}$ neighbors is thus given by

$g_{u}=\frac{1}{k_{u}}\sum_{v\in N(u)}g_{uv},$ (4)

being $N(u)$ the set of neighbors of $u$.

Evolution of the game is designed by means of the following update rule: at
each time step, all pairs of linked players engage in the PD game and their
payoffs $g_{u}$ are updated according to Eq. (4). Each player $u$ then
compares its performance with its neighborhood. Let us refer to the neighbor
of $u$ with the highest payoff as $w$, and as $x_{wu}$ its strategy adopted
against $u$. Then, $u$ updates its strategy $x_{uv}$ imitating its best
neighbor strategy $x_{wu}$ with probability

$p(\alpha,\Delta g_{u})=\frac{1}{1+e^{-\frac{\Delta g_{u}}{b\alpha}}},$ (5)

which is a Fermi function [69], i.e. a sigmoid function of the free parameter
$\alpha$ (playing the role of a temperature) and the payoff difference $\Delta
g_{u}=g_{w}-g_{u}$ (playing the role of an energy). Therefore, the update rule
of the strategy $x_{uv}(t)$ of $u$ at time $t$ is given by

$x_{uv}(t+1)=\left\\{\begin{array}[]{ll}x_{uv}(t)&\mbox{with probability
}1-p\\\ x_{wu}(t)&\mbox{with probability }p,\end{array}\right.$ (6)

that is, player $u$ imitates, at time $t+1$, the strategy that its best
neighbor $w$ played against it in the previous step with $k_{u}p$ of its
neighbors, and keeps the same strategy is actually playing with each one of
the rest.

One can define the node outcoming cooperation rate as

$\rho_{u}=\frac{1}{k_{u}}\sum_{v\in N(u)}x_{uv},$ (7)

and the incoming cooperation rate as

$\tau_{u}=\frac{1}{k_{u}}{\sum_{v\in N(u)}x_{vu}},$ (8)

and the macroscopic cooperation frequency of the whole network as

$\rho=\frac{1}{N}{\sum_{u=1}^{N}\rho_{u}}.$ (9)

Note that defection is very profitable when the benefit $b$ is close to $1$,
while when $b\to\infty$ the profitability of defection reduces. As a
consequence, in well-mixed populations, defection is the dominant trait
driving cooperation to extinction. However, it is known that, in structured
populations, cooperation can survive provided the appropriate conditions [70,
7].

In this context, we now explore how the transition from defection to
cooperation depends on the connectivity, as the dilemma weakness $b$
increases. To do this, let us consider two linked players $u$ and $v$ with
strategy vectors $\overline{x}_{uv}=(x_{uv},1-x_{uv})$ and
$\overline{x}_{vu}=(x_{vu},1-x_{vu})$. Using Eq. (2), one can write the payoff
of the first player from engaging the second as

$\overline{x}_{uv}P(b)\overline{x}_{vu}^{T}=x_{vu}-\frac{x_{uv}}{b}.$ (10)

Therefore, the total payoff accumulated by $u$ is

$g_{u}=\frac{1}{k_{u}}{\sum_{v\in
N(u)}\overline{x}_{uv}P(b)\overline{x}_{vu}^{T}}=\frac{1}{k_{u}}{\sum_{v\in
N(u)}\left(x_{vu}-\frac{x_{uv}}{b}\right)},$ (11)

or, equivalently, using definitions (7) and (8)

$g_{u}=\tau_{u}-\frac{\rho_{u}}{b}.$ (12)

Let us focus now on the probability of changing strategy, which becomes 50%
when the payoff difference $\Delta g_{u}=0$ in Eq. (5). Whenever at some point
during the iterated game a player $u$ reaches such a null payoff difference
with respect to its best performing neighbor $w$, the transition to keep the
actual strategies $x_{uv}$ or to change to $x_{wu}$ is critical in the sense
that it occurs for a particular value of the parameter $b$. This implies that,
evaluating payoffs of $u$ and $w$ using Eq. (12),

$\Delta g_{u}=\tau_{w}-\tau_{u}-\frac{1}{b*}(\rho_{w}-\rho_{u})=0,$ (13)

a critical value of $b^{*}$ is evaluated as

$b^{*}=\displaystyle\frac{\rho_{w}-\rho_{u}}{\tau_{w}-\tau_{u}}=\frac{\frac{1}{k_{w}}\sum_{z\in
N(w)}x_{wz}-\frac{1}{k_{u}}\sum_{v\in N(u)}x_{uv}}{\frac{1}{k_{w}}\sum_{z\in
N(w)}x_{zw}-\frac{1}{k_{u}}\sum_{v\in N(u)}x_{vu}.}$ (14)

From this relationship, one can formulate a simple conjecture based on the
connectivity pattern to locate the values of $b$ at which cooperation is
likely to transit to higher (or lower) frequencies: given a PD game on a graph
$G$, all phase transitions are located at points of the form

$b^{*}=\frac{p/k_{v}-q/k_{u}}{r/k_{v}-s/k_{u}}=\frac{pk_{u}-qk_{v}}{rk_{u}-sk_{v}},$
(15)

where $0\leq p,r\leq k_{v}$, $0\leq q,s\leq k_{u}$, are integers (such that
the denominator is non-zero).

Note that the set of $b^{*}$ values resulting from Eq. (15) tells us about the
possible locations of a phase transition. Equation (15) entirely depends on
the unique sequence of degrees present in the graph and, therefore, it alone
allows predicting the transitions without actually playing the game. However,
the multiplicity of each $b^{*}$ does depend on the particular degree
distribution, and one can extract important information about the transitions
from this degeneracy, as shown in the next section.

While in the game framework described so far, a player can perform different
strategies with its neighbors, in the more classic scenario, players adopt a
single role in each run. We will refers to this updating rule as node
strategy, in opposition to the link strategy defined by Eq. (6). It is
straightforward to obtain the equivalent equation to that of Eq. (15) in this
latter case:

$b^{*}=\frac{1}{r/k_{v}-s/k_{u}}=\frac{k_{u}k_{v}}{rk_{v}-sk_{u}},$ (16)

where $0\leq r\leq k_{u}$ and $0\leq s\leq k_{v}$ are integers (such that
denominator is non-zero) and Eq. (16) indicates that the possible locations
for a phase transition in a node strategy situation is reduced with respect to
the link strategy choice.

Notice that the two previous conjectures for link ( mixed) and node
strategies, Eqs. (15) and (16) respectively, depend on the particular form of
the payoff matrix. Therefore, the relationship between game parameters and
topology has to be recomputed if we change the dilemma, as the critical points
$\Delta g_{u}=0$ depend on them. For example, in order to illustrate this with
a case extensively investigated in the literature, we chose the dilemma’s
matrix with no cost

$P^{\prime}(b)=\begin{pmatrix}1&0\\\ b&0\\\ \end{pmatrix}.$ (17)

as in Refs.[7, 71, 72]. Note that parameter $b$ in Eq. (17) now plays the role
of a temptation to defect. In this case, the corresponding universal dilemma
strengths [59] are $D^{\prime}_{g}=b$ and $D^{\prime}_{r}=0$, therefore,
parameter $b$ directly accounts for the strength of the dilemma.

Using a similar procedure as the one leading to Eq. (16) one can obtain an
equivalent conjecture to locate the phase transitions for the dilemma
$P^{\prime}$ above,

$b^{*}=\frac{n}{m}$ (18)

where $n$ and $m$ are integers such that $1\leq m\leq n\leq\max(k)$, and
$\max(k)$ is the maximum degree present in the network.

##  3 Results

Figure 1:  Dependence of the conjecture on the set of different degrees
present in the network. Average cooperation frequency $\rho$ as a function of
the dilemma weakness $b$ for random regular graphs of size $N=600$ and
different degree distributions: (a) $k=\\{5\\}$, (b) $k=\\{3,5\\}$
($N_{3}=N_{5}=300$), (c) $k=\\{2,3,5\\}$ ($N_{2}=N_{3}=N_{5}=200$). The
vertical grey dotted lines correspond to the $b^{*}$ values computed with Eq.
(15). The curves are averaged over 100 network realizations. Other parameter
values: $\alpha=0.05$, $T=10^{4}$.

We test our conjecture by running an extensive set of simulations of PD
competition in a wide range of structural and dynamical conditions. Unless
otherwise indicated, all the simulations were carried out for a population of
$N=600$ players using link strategy, with initially equal densities of
cooperation and defection randomly distributed along with the links. Each
simulation evolves over 10,000 time steps, and cooperation frequencies are
averaged in the last 2,500 time steps. At each time step, all players update
their strategies synchronously [1, 73]. Average cooperation frequency is
monitored as a function of $b$, and each point is an ensemble average over 100
network realisations.

Figure 2:  Dependence of the conjecture on the heterogeneity of the degree
distribution and multiplicity of the critical game parameter. (Top panels)
Average cooperation frequency $\rho$ as a function of the dilemma weakness $b$
for ER (left) and SF (right) graphs of size $N=600$ and $\langle k\rangle=4$.
Each curve is the average of 100 simulations of the same type of network.
(Bottom panels) Probability density function of the $b^{*}$ resulting from Eq.
(15) for the ER (left) and SF (right) degree distributions. Vertical dotted
lines in all panels correspond to the local maxima whose prominence is above
the threshold marked by the red horizontal line ($20\%$ of the maximum value).

As the conjecture implies that the number of possible transitions depends on
the set of degrees present in the network, we initiate our exploration by
showing the results of simulations on random regular (RR) graphs composed of
groups of equal size $N_{k}$ whose nodes have all the same connectivity $k$.
These graphs are constructed using the configuration model [74], which returns
a random graph with no degree correlations consistent with a given degree
distribution. Specifically, we explore in Fig. 1 the evolution of the
cooperation frequency in RR graphs for the following degree distributions: a
unique group of $N_{k=5}=N$ nodes of degree $k={5}$ (upper panel); two equally
sized groups $N_{3}=N_{5}=N/2$ of agents with degrees $k=\\{3,5\\}$ (middle
panel); and three equally sized groups $N_{2}=N_{3}=N_{5}=N/3$ with degrees
$k=\\{2,3,5\\}$ respectively (bottom panel). In addition, we have included
dotted vertical lines marking the values of $b^{*}$ predicted by Eq. (15) for
the chosen degree sequences. Note that the maximum degree is the same in the
three cases but, as smaller degrees are included, the cardinality of $b^{*}$
increases adding new locations for a potential phase transition. We observe
that all significant changes in the cooperation frequency always occur at one
of the possible $b^{*}$ values predicted by Eq. (15).

Figure 3:  Dependence of the conjecture on the network assortativity. (Top
row) Average cooperation frequency $\rho$ versus the dilemma weakness $b$ of
random regular graphs with degrees $k=\\{2,3\\}$ for (left) different
population size and (right) different levels of assortativity. (Bottom row)
Frequency distribution of the corresponding $b^{*}$ values calculated
considering the actual degree pairings in the networks.

In order to inspect the role of the degree distribution in more detail, we
carried out simulations with random Erdös-Renyi [75] and scale-free Barabási-
Albert [76] networks, with $N=600$ and average degree $\langle k\rangle=4$. In
the upper panels of Fig. 2 we observe that ER networks support higher values
of the cooperation frequency and that in both cases, cooperation increases in
a stepwise manner, with more pronounced plateaux in the ER case. For such more
complex networks, the set of possible $b^{*}$ values predicted by Eq.(15) is
very numerous, and one sees that not in all of them does a significant change
in the cooperation frequency appears. To explain this, we consider not only
the sequence of unique degrees, but the degree distribution. By knowing how
often each degree $k$ appears in the network, we can calculate the
multiplicity of each $b^{*}$ value, resulting in an excellent indicator of how
probable is that a certain value $b^{*}$ determines a sharp variation in the
cooperation frequency. In the bottom panels of Fig. 2 we plot the histograms
of the average multiplicity of the $b^{*}$ values in Eq. (15) in the degree
distributions of the ER (left) and SF (right) networks. The vertical dotted
lines mark those local maxima whose prominence (height difference with respect
to the local surrounding baseline) is above a given threshold (red horizontal
line). Remarkably, the steps transitions observed in the cooperation frequency
occur precisely at these more frequent $b^{*}$, as shown by the same vertical
dotted lines in the top panels.

To further investigate the role played by the second order structural details
we study the impact of varying the correlation degree properties of the
network while keeping its degree sequence. Figure 3 (top left panel) reports
how cooperation evolves in RR networks which include only nodes with degrees
$k=\\{2,3\\}$, but controlling the population size $N_{k}$: equal size
($N_{2}=N_{3}$) and two asymmetric cases where one population is much larger
than the other ($N_{3}=5N_{2}$ or $N_{2}=5N_{3}$). One clearly observes that
the relevant jumps in the cooperation occur at exactly the same values of $b$,
namely, $b=2,3,6$, matching the predicted $b^{*}$ values highlighted by the
vertical dotted lines. However, the steepness of the transitions is certainly
not the same due to a difference in the $b^{*}$ multiplicity caused by
different degree frequencies.

A more subtle but showing case is presented in Figure 3 (top right panel).
Here the degree distribution is preserved in all cases and $N_{2}=N_{3}=300$.
Therefore, the theoretical $b^{*}$ multiplicity coming from the degree
distribution is the same. However, in this example we control the level of
degree-degree correlation measured by the associated Pearson coefficient $r$
[77] (see legend in the right panel). Starting from uncorrelated networks
($r=0$), we perform the rewiring procedure of Ref. [78] that preserves the
degree distribution until the desired level of degree-degree correlation is
achieved [79, 80]. This way, one can produce graphs with very high positive
correlation, in which nodes with the same degree tend to be linked among them
($r\sim 1$), and networks whose nodes have $k=2$ are more likely to be linked
to nodes with $k=3$ ($r\sim-0.7$). The clear differences in the transitions
(in particular for positively correlated networks) indicates that, since the
calculation of $b^{*}$ involves all possible pairs of degrees, in correlated
networks the real multiplicity could be shifted from the neutral one
determined only by the degree distribution. Indeed, the bottom panels in Fig.
3 show the histograms for the actual multiplicity of $b^{*}>k_{\rm min}$
considering the real pairings in the adjacency matrices, revealing that the
height of the bars differs in each case, and they provide a closer predictor
of a phase transition in the game.

Figure 4:  Assessment of the impact of the temperature on the conjecture’s
predictive performance. Cooperation dependence on the temperature $\alpha$ for
RR graphs ($k=\\{3,5\\}$ and $N_{3}=N_{5}$). Panel (a) shows the cooperation
$\rho$ as a function of $b$ for several values of $\alpha$. Vertical dotted
lines are the corresponding $b^{*}$ values predicted by Eq. (15) for the RR
network. Panel (b) shows the evolution as a function of $\rho$ for values of
$b$ chosen close to a transition in (a). The network size is $N=600$ and each
curve is an average over 100 network realizations. Figure 5:  Extension of the
conjecture to the prisoner’s dilemma game on networks using node strategy.
Average cooperation frequency $\rho$ as a function of the dilemma weakness $b$
in a PD game using node strategy on different networks of size $N=600$ and
$\langle k\rangle=4$: RR with $k=\\{3,5\\}$ ($N_{3}=N_{5}$), ER and SF. Curves
are averaged over 100 network realizations. Vertical lines correspond to the
most frequent $b^{*}$ values given by Eq. (16) after substituting the degree
sequences of the RR (blue dashed), ER (red dotted) and SF (orange dotted)
networks. Other parameter values: $\alpha=0.05$.

Notice that our conjecture states that if there is a sharp change in the
cooperation frequency, the corresponding $b$ is one of the $b^{*}$ values, but
the reverse is not always true. This unidirectional correspondence can be
verified, for instance, in Fig. 1, where not all the predicted $b^{*}$ values
(indicated by the dotted vertical lines) lead to an abrupt change in the
cooperation. As the conjecture involves a critical phenomenon, we also have to
carefully inspect the role of the stochastic part of the dynamics, here
represented by the temperature parameter $\alpha$ which controls the shape of
the probability function for the strategy update. To this end, in Fig. 4 we
show the evolution of cooperation in RR graphs with $k=\\{3,5\\}$
($N_{3}=N_{5}$ as in Fig. 1(b)) for several values of the temperature. In
addition, we also include the corresponding $b^{*}$ values as vertical dotted
lines. One easily sees substantial changes in cooperation always coinciding
with one of the vertical lines but, depending on the temperature, different
$b^{*}$ are selected, and also differences in the height of the transition are
observed. For instance, while lower temperatures support cooperation changes
around $b^{*}=7.5$, for $\alpha=0.1$, cooperation does not evolve. Temperature
dependence of the cooperation is shown on the right panel for values of $b$
close to the major transitions observed in the left panel. The results suggest
the existence of an optimal temperature that maximises cooperation, with a
peak that shifts to lower temperatures as $b$ increases.

The overall scenario just described in the case of link strategy remains
qualitatively unchanged for the simpler case of node strategy, considered as a
particular case in Eq. (16). Figure 5 shows the cooperation curves for node
strategy in RR graphs with $k=\\{3,5\\}$ and ER and SF graphs with the same
average degree $\langle k\rangle=4$. In general, as compared with Fig. 1(b)
and Figs. 2(a,b), one observes that cooperation frequencies reach higher
values with node strategy games. Although Eq. (16) predicts a smaller set of
$b^{*}$ than Eq.(15) for the same degree sequence, it is still able to capture
the major transitions as shown by the vertical lines corresponding to the
$b^{*}$ with higher prominence in their frequency distribution. For example,
if we closely compare the transition to cooperation in the RR graphs using
link and node strategies, Fig. 1(b) and Fig.5 (blue line), major changes take
place at $b=3.0,3.75$ and $5$, as predicted by both schemes.

For the sake of comparison with other classic studies [7, 18], we also
numerically check the conjecture obtained in Eq. (18) for the unpunished case
$P^{\prime}$ introduced at the end of Sec. 2. Furthermore, in this case the
system evolution is carried out by implementing Monte Carlo simulations such
that whenever the payoffs of two engaged $u$ and $v$ randomly selected players
verify $g_{v}>g_{u}$, then player $u$ imitates $v$’s strategy with a
probability $p=\frac{D}{2\Gamma}$ for $D<2\Gamma$ and $p=1$ otherwise, being
$D=(g_{v}-g_{u})/(b\max(k_{u},k_{v}))$ a normalized payoff difference. Notice
that the slope $\Gamma$ is playing a similar functional role than $\alpha$ in
the sigmoidal function of Eq. (5) and that $\Gamma=0.5$ corresponds to the
case reported in Refs. [81, 18].

Figure 6:  Extension of the conjecture to the prisoner’s dilemma game with no
cost and non synchronous strategy update of randomly selected agents. Average
cooperation frequency $\rho$ as a function of the temptation to defect (or
dilemma strength) $b$ for ER (left) and SF (right) graphs of size $N=500$ and
$\langle k\rangle=4$. Each curve is the average of 100 simulations of the same
type of network implementing a Monte Carlo method following a transition
probability given by $p=\frac{D}{2\Gamma}$ for $0\leq D\leq 2\Gamma$ and $p=1$
for $D>2\Gamma$. Equilibrium frequencies are obtained by averaging in the last
2,500 generations before the transient time of 10,000 generations. Vertical
dotted lines in both panels correspond to those $b^{*}$ values from Eq. (18)
whose pdf is above 10% and 5% for ER and SF networks respectively.

Figure 6 shows how the cooperation evolves as a function of the game parameter
$b$ for the dilemma with no cost for both ER and SF networks and two values of
the parameter $\Gamma$ quantifying the uncertainty of strategy adoptions. We
observe that cooperation evolves in a continuous way for $\Gamma=0.5$
(transition probability is directly the normalised payoff difference as used
in Ref. [7]) and, phase transitions from pure cooperation to a mixed state
(and from a mixed state to complete defection for ER topologies) are perfectly
predicted by those $b^{*}$ values (marked as dotted vertical lines) in Eq.
(18) whose probability density function is above a given threshold. One sees
that reducing the value of $\Gamma$ (i.e. increasing the likelihood of a
player to change its strategy even for very low payoff differences) breaks the
continuous evolution of the cooperation into a series of plateaux at different
mixed levels of defection and cooperation.

##  4 Discussion

Cooperation is much more widespread in nature than the Darwinian premise of
‘only the fittest survive’ might suggest. Why and under which conditions
cooperation thrives is, therefore, an evergreen subject across the social and
natural sciences [82, 83, 84, 85, 86, 87]. Evolutionary game theory is
traditionally used to formalise the problem mathematically with social
dilemmas [88, 89, 90], and networks are commonly used as the backbone for the
simulation of these games. We have here shown that the structure of the former
in terms of the unique sequence of degrees predict the outcome in the latter,
particularly in terms of the game parameters at which significant shifts in
the level of cooperation can be expected. In particular, we have proposed a
conjecture for cooperation transitions in arbitrary networks, including phase
transitions from absorbing to mixed strategy phases.

Results based on the evolutionary prisoner’s dilemma game on random and scale-
free networks of various sizes demonstrate the effectiveness of the
conjecture. Indeed, since the conjecture is based on the different degrees
present in the network, the network size does not modify the set of game
parameters critical to produce shifts in the cooperation frequency. We have
reported good agreements between the predictions and the results obtained with
either a synchronous and asynchronous update of strategies, thus providing a
viable alternative to computationally expensive and often time-consuming
simulations. We have shown that beyond the main qualitative information
provided by the unique sequence of degrees, the degree distribution and the
possible second order structural correlations provide further quantitative
details about the transitions. Although we have explored the predictive power
of the conjecture only for the prisioner’s dilemma, it can be generalized to
any social dilemma as the conjecture just operates on the payoff matrix and
the critical game parameters are extracted only from the condition that the
payoff difference between two agents is zero. Despite the success, we also
note the requirements for the conjecture to work best. Namely, the presence of
high stochasticity levels in the probability of strategy update may compromise
the emergence of sharp transitions, which can explain the fact that these have
not been previously observed.

In this light, the conjecture holds the promise of the optimisation of large-
scale simulations of evolutionary processes on complex networks, not just
related to cooperation as a particular example of prosocial behaviour, but
more generally for the broader class of moral behaviour where similar models
are often employed [91]. By providing a fast and straightforward way to
estimate the outcome of evolutionary processes on arbitrary networks and a
wide class of games without the need for actual simulations, the conjecture
should prove helpful in a wide range of research fields, further strengthening
the role of physics and network science within them.

I.S.N. and I.L. acknowledge support from the Ministerio de Economía, Industria
y Competitividad of Spain under project FIS2017-84151-P and the Ministerio de
Ciencia e Innovación under project PID2020-113737GB-I00. They also acknowledge
the computational resources and assistance provided by CRESCO, the
supercomputing center of ENEA in Portici, Italy. M.P. was supported by the
Slovenian Research Agency (Grant Nos. P1-0403 and J1-2457).

## References

## References

  * [1] Nowak M A and May R M 1992 Nature 359 826–829
  * [2] Nash J 1950 Proc. Natl. Acad. Sci. U.S.A. 36 48–49
  * [3] Nash J 1951 Ann. Math. 54 286–295
  * [4] Szabó G and Fáth G 2007 Phys. Rep. 446 97–216
  * [5] Perc M and Szolnoki A 2010 BioSystems 99 109–125
  * [6] Wang Z, Wang L, Szolnoki A and Perc M 2015 Eur. Phys. J. B 88 124
  * [7] Santos F C and Pacheco J M 2005 Phys. Rev. Lett. 95 098104
  * [8] Nowak M A 2006 Science 314 1560–1563
  * [9] Albert R and Barabási A L 2002 Rev. Mod. Phys. 74 47–97
  * [10] Newman M E J 2003 SIAM Review 45 167–256
  * [11] Boccaletti S, Latora V, Moreno Y, Chavez M and Hwang D 2006 Phys. Rep. 424 175–308
  * [12] Santos F C, Pacheco J M and Lenaerts T 2006 Proc. Natl. Acad. Sci. USA 103 3490–3494
  * [13] Gómez-Gardeñes J, Campillo M, Floría L M and Moreno Y 2007 Phys. Rev. Lett. 98 108103
  * [14] Rong Z, Li X and Wang X 2007 Phys. Rev. E 76 027101
  * [15] Poncela J, Gómez-Gardeñes J, Floría L M and Moreno Y 2007 New J. Phys. 9 184
  * [16] Masuda N 2007 Proc. R. Soc. B 274 1815–1821
  * [17] Tomassini M, Luthi L and Pestelacci E 2007 Int. J. Mod. Phys. C 18 1173–1185
  * [18] Szolnoki A, Perc M and Danku Z 2008 Physica A 387 2075–2082
  * [19] Assenza S, Gómez-Gardeñes J and Latora V 2008 Phys. Rev. E 78 017101
  * [20] Santos F C, Santos M D and Pacheco J M 2008 Nature 454 213–216
  * [21] Poncela J, Gómez-Gardeñes J and Moreno Y 2011 Phys. Rev. E 83 057101
  * [22] Tanimoto J, Brede M and Yamauchi A 2012 Phys. Rev. E 85 032101
  * [23] Amaral M A and Javarone M A 2018 Phys. Rev. E 97 042305
  * [24] Abramson G and Kuperman M 2001 Phys. Rev. E 63 030901(R)
  * [25] Kim B J, Trusina A, Holme P, Minnhagen P, Chung J S and Choi M Y 2002 Phys. Rev. E 66 021907
  * [26] Masuda N and Aihara K 2003 Phys. Lett. A 313 55–61
  * [27] Santos F C, Rodrigues J F and Pacheco J M 2005 Phys. Rev. E 72 056128
  * [28] Tomassini M, Luthi L and Giacobini M 2006 Phys. Rev. E 73 016132
  * [29] Fu F, Liu L H and Wang L 2007 Eur. Phys. J. B 56 367–372
  * [30] Vukov J, Szabó G and Szolnoki A 2008 Phys. Rev. E 77 026109
  * [31] Vukov J and Szabó G 2005 Phys. Rev. E 71 036133
  * [32] Lee S, Holme P and Wu Z X 2011 Phys. Rev. Lett. 106 028702
  * [33] Zimmermann M G, Eguíluz V M and San Miguel M 2004 Phys. Rev. E 69 065102(R)
  * [34] Pacheco J M, Traulsen A and Nowak M A 2006 Phys. Rev. Lett. 97 258103
  * [35] Fu F, Chen X, Liu L and Wang L 2007 Physica A 383 651–659
  * [36] Szolnoki A, Perc M and Danku Z 2008 EPL 84 50007
  * [37] Fu F, Wu T and Wang L 2009 Phys. Rev. E 79 036101
  * [38] Szolnoki A and Perc M 2009 New J. Phys. 11 093033
  * [39] Holme P, Trusina A, Kim B J and Minnhagen P 2003 Phys. Rev. E 68 030901
  * [40] Guimerà R, Danon L, Díaz-Guilera A, Giralt F and Arenas A 2003 Phys. Rev. E 68 065103
  * [41] Lozano S, Arenas A and Sánchez A 2008 PLoS ONE 3 e1892
  * [42] Boccaletti S, Bianconi G, Criado R, del Genio C, Gómez-Gardeñes J, Romance M, Sendiña-Nadal I, Wang Z and Zanin M 2014 Phys. Rep. 544 1–122
  * [43] Kivelä M, Arenas A, Barthelemy M, Gleeson J P, Moreno Y and Porter M A 2014 J. Complex Netw. 2 203–271
  * [44] Gómez-Gardeñes J, Gracia-Lázaro C, Floría L M and Moreno Y 2012 Phys. Rev. E 86 056113
  * [45] Jiang L L and Perc M 2013 Sci. Rep. 3 2483
  * [46] Wang Z, Wang L and Perc M 2014 Phys. Rev. E 89 052813
  * [47] Battiston F, Perc M and Latora V 2017 New J. Phys. 19 073017
  * [48] Shen C, Chu C, Shi L, Jusup M, Perc M and Wang Z 2018 EPL (Europhysics Letters) 124 48003
  * [49] Fotouhi B, Momeni N, Allen B and Nowak M A 2018 Nat. Human Behav. 2 492
  * [50] Shi L, Shen C, Geng Y, Chu C, Meng H, Perc M, Boccaletti S and Wang Z 2019 Nonlinear Dyn. 96 49–56
  * [51] Duh M, Gosak M, Slavinec M and Perc M 2019 New J. Phys. 21 123016
  * [52] Sun C, Luo C and Li J 2020 Chaos, Solitons & Fractals 135 109769
  * [53] Wang Z, Szolnoki A and Perc M 2013 Sci. Rep. 3 1183
  * [54] Battiston F, Cencetti G, Iacopini I, Latora V, Lucas M, Patania A, Young J G and Petri G 2020 Phys. Rep. 874 1–92
  * [55] Burgio G, Matamalas J T, Gómez S and Arenas A 2020 Entropy 22 744
  * [56] Alvarez-Rodriguez U, Battiston F, de Arruda G F, Moreno Y, Perc M and Latora V 2021 Nat. Hum. Behav. 5 586–595
  * [57] Guo H, Jia D, Sendiña-Nadal I, Zhang M, Wang Z, Li X, Alfaro-Bittner K, Moreno Y and Boccaletti S 2021 Chaos, Solitons & Fractals 150 111103 ISSN 0960-0779
  * [58] Tanimoto J and Sagara H 2007 Biosystems 90 105–114 ISSN 0303-2647
  * [59] Wang Z, Kokubo S, Jusup M and Tanimoto J 2015 Physics of Life Reviews 14 1–30 ISSN 1571-0645
  * [60] Ito H and Tanimoto J 2018 Royal Society Open Science 5 181085
  * [61] Ohtsuki H, Hauert C, Lieberman E and Nowak M A 2006 Nature 441 502–505
  * [62] Allen B, Lippner G, Chen Y T, Fotouhi B, Nowak M A and Yau S T 2017 Nature 544 227–230
  * [63] Javarone M A, Antonioni A and Caravelli F 2016 EPL 114 38001
  * [64] Javarone M A 2016 Eur. Phys. J. B 89 42
  * [65] Rapoport A 1989 Prisoner’s dilemma Game Theory (Springer) pp 199–204
  * [66] Sendiña Nadal I, Leyva I, Perc M, Papo D, Jusup M, Wang Z, Almendral J A, Manshour P and Boccaletti S 2020 Phys. Rev. Research 2(4) 043168
  * [67] Kokubo S, Wang Z and Tanimoto J 2015 Appl. Math. Comput. 259 552–568
  * [68] Tanimoto J 2017 Appl. Math. Comput. 304 20–27
  * [69] Fermi E 1926 Rendiconti Lincei 3 145–149
  * [70] Ohtsuki H, Hauert C, Lieberman E and Nowak M A 2006 Nature 441 502–505
  * [71] Nowak M A and Sigmund K 1992 Nature 355 250–253
  * [72] Perc M 2009 New J. Phys. 11 033027
  * [73] Huberman B and Glance N 1993 Proc. Natl. Acad. Sci. USA 90 7716–7718
  * [74] Molloy M and Reed B 1995 Random Structures & Algorithms 6 161–180
  * [75] Erdös P and Rényi A 1960 Publ. Math. Inst. Hung. Acad. Sci 5 17–60
  * [76] Barabási A L and Albert R 1999 science 286 509–512
  * [77] Newman M E J 2003 Phys. Rev. E 67(2) 026126
  * [78] Sendiña Nadal I, Leyva I, Navas A, Villacorta-Atienza J A, Almendral J A, Wang Z and Boccaletti S 2015 Phys. Rev. E 91(3) 032811
  * [79] Maslov S and K S 2002 Science 296 910–3
  * [80] Xulvi-Brunet R and Sokolov I M 2004 Phys. Rev. E 70(6) 066102
  * [81] Santos F C and Pacheco J M 2005 Phys. Rev. Lett. 95(9) 098104
  * [82] Axelrod R 1984 The Evolution of Cooperation (New York: Basic Books)
  * [83] Henrich N and Henrich J P 2007 Why Humans Cooperate: A Cultural and Evolutionary Explanation (Oxford, U.K.: Oxford University Press)
  * [84] Nowak M A and Highfield R 2011 SuperCooperators: Altruism, Evolution, and Why We Need Each Other to Succeed (New York: Free Press)
  * [85] Rand D G and Nowak M A 2013 Trends in Cognitive Sciences 17 413–425
  * [86] Kraft-Todd G, Yoeli E, Bhanot S and Rand D 2015 Current Opinion in Behavioral Sciences 3 96–101
  * [87] Perc M, Jordan J J, Rand D G, Wang Z, Boccaletti S and Szolnoki A 2017 Phys. Rep. 687 1–51
  * [88] Weibull J W 1995 Evolutionary Game Theory (Cambridge, MA: MIT Press)
  * [89] Hofbauer J and Sigmund K 1998 Evolutionary Games and Population Dynamics (Cambridge, U.K.: Cambridge University Press)
  * [90] Nowak M A 2006 Evolutionary Dynamics (Cambridge, MA: Harvard University Press)
  * [91] Capraro V and Perc M 2021 J. R. Soc. Interface 18 20200880