# A New Limit on Planck Scale Lorentz Violation from Gamma Ray Burst
Polarization

Floyd W. Stecker  Astrophysics Science Division, NASA Goddard Space Flight
Center, Greenbelt, MD 20771, USA<EMAIL_ADDRESS>
###### Abstract

Constraints on possible Lorentz invariance violation (LIV) to first order in
$E/M_{\rm Planck}$ for photons in the framework of effective field theory
(EFT) are discussed, taking cosmological factors into account. Then, using the
reported detection of polarized soft $\gamma$-ray emission from the
$\gamma$-ray burst GRB041219a that is indicative of an absence of vacuum
birefringence, together with a very recent improved method for estimating the
redshift of the burst, we derive constraints on the dimension 5 Lorentz
violating modification to the Lagrangian of an effective local QFT for QED.
Our new constraints are more than five orders of magnitude better than recent
constraints from observations of the Crab Nebula. We obtain the upper limit on
the Lorentz violating dimension 5 EFT parameter $|\xi|$ of $2.4\times
10^{-15}$, corresponding to a constraint on the dimension 5 standard model
extension parameter, $k^{(5)}_{(V)00}\leq 4.2\times 10^{-34}$ GeV-1.

###### keywords:

Lorentz invariance; quantum gravity; gamma-rays; gamma-ray bursts

††journal: Astroparticle Physics

##  1 Introduction

Because of the problems associated with merging relativity with quantum
theory, it has long been felt that relativity will have to be modified in some
way in order to construct a quantum theory of gravitation. Since the Lorentz
group is unbounded at the high boost (or high energy) end, in principle it may
be subject to modifications in the high boost limit [1, 2]. There is also a
fundamental relationship between the Lorentz transformation group and the
assumption that space-time is scale-free, since there is no fundamental length
scale associated with the Lorentz group. However, as noted by Planck [3],
there is a potentially fundamental scale associated with gravity, viz., the
Planck scale $\lambda_{Pl}=\sqrt{G\hbar/c^{3}}\sim 10^{-35}$ m, corresponding
to an energy (mass) scale of $M_{Pl}=\hbar c/\lambda_{Pl}\sim 10^{19}$ GeV.

In recent years, there has been much interest in testing Lorentz invariance
violating terms that are of first order in $E/M_{Pl}$, since such terms vanish
at very low energy and are amenable to testing at higher energies. In
particular, tests using high energy astrophysics data have proved useful in
providing constraints on Lorentz invariance violation (LIV) (e.g., see reviews
in Refs. [4] and [5]).

##  2 Vacuum Birefringence

Important fundamental constraints on LIV come from searches for the vacuum
birefringence effect predicted within the framework of the effective field
theory (EFT) analysis of [6]. (See also Ref. [7]). Within this framework,
applying the Bianchi identities to the leading order Maxwell equations in
vacua, a mass dimension 5 operator term is derived of the form

${\Delta\cal{L}}_{\gamma}={{\xi}\over{M_{Pl}}}{n^{a}F_{ad}n\cdot\partial(n_{b}\tilde{F}^{bd})}.$
(1)

It is shown in Ref. [6] that the expression given in Equation (1) is the only
dimension 5 modification of the free photon Lagrangian that preserves both
rotational symmetry and gauge invariance. This leads to a modification in the
dispersion relation proportional to $\xi(\omega/M_{Pl})=\xi(E/M_{Pl})$
111adopting the conventions $\hbar=1$ and the low energy speed of light $c=1$.
with the new dispersion relation given by

$\omega^{2}~{}=~{}k^{2}\pm\xi\,k^{3}/M_{Pl}.\\\ $ (2)

Some models of quantized space-time suggest $\xi$ should be $\cal{O}$(1),
(see, e.g., Ref. [8]). The sign in the photon dispersion relation corresponds
to the helicity, i.e., right or left circular polarization. Equation (2)
indicates that photons of opposite circular polarization have different phase
velocities and therefore travel with different speeds. The effect on photons
from a distant linearly polarized source can be constructed by decomposing the
linear polarization into left and right circularly polarized states. It is
then apparent that this leads to a rotation of the linear polarization
direction through an angle

$\theta(t)=\left[\omega_{+}(k)-\omega_{-}(k)\right]t_{P}/2~{}\simeq~{}\xi
k^{2}t_{P}/2M_{Pl}$ (3)

for a plane wave with wave-vector $k$, where $\xi k/M_{Pl}\ll 1$ and where
$t_{P}$ is the propagation time.

Observations of polarized radiation from distant sources can thus be used to
place an upper bound on $\xi$. The vacuum birefringence constraint arises from
the fact that if the angle of polarization rotation (3) were to differ by more
than $\pi/2$ over the energy range covered by the observation the
instantaneous polarization at the detector would fluctuate sufficiently for
the net polarization of the signal to be suppressed well below any observed
value. The difference in rotation angles for wave-vectors $k_{1}$ and $k_{2}$
is

$\Delta\theta=\xi(k_{2}^{2}-k_{1}^{2})L_{P}/2M_{Pl},$ (4)

where we have replaced the propagation time $t_{P}$ by the propagation
distance $L_{P}$ from the source to the detector.

If polarization is detected from a source at redshift $z$, this yields the
constraint

$|\xi|<{{\pi
M_{Pl}}\over{\int\displaylimits_{0}^{z}dz^{\prime}{[k_{2}(z^{\prime})^{2}-k_{1}(z^{\prime})^{2}]|dL_{P}(z^{\prime})/dz^{\prime}|}}}$
(5)

where $k_{1,2}(z^{\prime})=(1+z^{\prime})\cdot k_{1,2}(z^{\prime}=0)$ and

$\Bigl{|}{dL_{P}\over{dz^{\prime}}}\Bigr{|}={{c}\over{H_{0}}}{{1}\over{(1+z^{\prime})\sqrt{\Omega_{\Lambda}+(1+z^{\prime})^{3}\Omega_{m}}}}.$
(6)

Defining

${\cal{D}}={{c}\over{H_{0}}}\int\displaylimits_{0}^{z}dz^{\prime}{{(1+z^{\prime})}\over{\sqrt{\Omega_{\Lambda}+(1+z^{\prime})^{3}\Omega_{m}}}}$
(7)

it follows from equations (5)-(7) and the definitions of $k_{1,2}(z^{\prime})$
that

$|\xi|<{{\pi M_{Pl}}\over{{\cal{D}}(k_{2}^{2}-k_{1}^{2})}},$ (8)

with the standard cosmological values [9] of $\Omega_{m}=0.27$,
$\Omega_{\Lambda}=0.73$, and $H_{0}$ = 71 km s-1 Mpc-1 (1 Mpc = $3.09\times
10^{22}$ m). Figure 1 shows the function ${\cal{D}}(z)$ as defined in Equation
(7).

Figure 1: A linear plot of the integral $\cal{D}$ as defined in Equation (7),
given as a function of redshift, $z$.

##  3 Previous Constraints

A previous bound of $|\xi|\lesssim 2\times 10^{-4}$, was obtained by Gleiser
and Kozameh [10] using the observed 10% polarization of ultraviolet light from
a galaxy at distance of around 300 Mpc. Fan et al. used the observation of
polarized UV and optical radiation at several wavelengths from the
$\gamma$-ray bursts (GRBs) GRB020813 at a redshift $z=1.3$ and GRB021004
$z=2.3$ to get a constraint of $|\xi|\lesssim 2\times 10^{-7}$ [11]. Jacobson
et al. [12] used a report of polarized $\gamma$-rays observed [13] in the
prompt emission from the $\gamma$-ray burst GRB021206 in the energy range 0.15
to 2 MeV using the RHESSI detector [14] to place strong limits on $\xi$.
However, this claimed polarization detection has been refuted [15, 16].

Kostelecký and Mewes [17] have shown that the EFT model parameter $\xi$ can be
related to the model independent isotropic dimension 5 standard model
extension (SME) parameter $k^{(5)}_{(V)00}$. They derive the relation

$k^{(5)}_{(V)00}=3\sqrt{4\pi}\xi/5M_{Pl},$ (9)

which we use in this paper. Their upper limit of $1\times 10^{-32}$ GeV-1,
obtained by assuming a lower limit on the redshift of these bursts of $z=0.1$,
then corresponds to the constraint $\xi<6\times 10^{-14}$.222Ref. [18] gives a
table of similar limits on $k^{(5)}_{(V)00}$ with citations.

More recently, Maccione et al. have derived a constraint of $|\xi|\lesssim
9\times 10^{-10}$ using observations of polarized hard X-rays from the Crab
Nebula detected by the INTEGRAL satellite [19].

It is clear from Equation (5) that the larger the distance of the polarized
source, and the larger the energy of the photons from the source, the greater
the sensitivity to small values of $\xi$. In that respect, the ideal source to
study would be polarized X-rays or $\gamma$-rays from a GRB with a known
redshift at a deep cosmological distance [12].

##  4 A New Treatment

Unfortunately, despite the many GRBs that have been detected and have known
host galaxy spectral redshifts, none of these bursts have measured
$\gamma$-ray polarization. However, in this paper we take a new approach,
deriving an estimated redshift for GRB041219a. This is a GRB with reported
polarization but no spectral redshift measurement.

Polarization at a level of 63(+31,-30)% to 96(+39,-40)% in the soft
$\gamma$-ray energy range has been detected by analyzing data from the
spectrometer on INTEGRAL for GRB041219a in the 100 to 350 keV energy range
[20]. It should be noted that that a systematic effect that might mimic
polarization in the analysis could not definitively be excluded. This GRB does
not have an associated host galaxy spectral redshift.

Useful relations have been recently obtained where known spectral redshifts of
GRBs are statistically correlated with various observational parameters of the
bursts such as luminosity, the Band function [21] parameter $E_{peak}$, rise
time, lag time and variability of a burst (Ref. [22] and references therein).
A detailed treatment of these correlations is given in Ref. [22]. By deriving
updated luminosity correlations for a very large number of GRBs, they find the
tightest correlation is the luminosity-$E_{peak}$ correlation. Using the
relation given in Ref. [22],

$\log L=51.75+1.35\log[(1+z)E_{peak}/300{\rm keV}]$ (10)

and the iterative method described in Ref. [23], and taking $E_{peak}$ = 170
keV and a peak fluance of $5.7\times 10^{-4}$ erg cm-2 [20], we derive a value
for $z$ for GRB041219a of $0.23\pm 0.03$. Taking a lower limit of 0.2 for the
redshift and taking $k_{2}$ = 350 keV/c and $k_{1}$ = 100 keV/c in Equation
(5), we find a new, most accurate cosmological constraint on $|\xi|$ of

$|\xi|\leq 2.4\times 10^{-15},$ (11)

almost five orders of magnitude better than the previous best solid limit
derived using polarimetric observations of the Crab Nebula in the hard X-ray
energy range [19].

From equation 9, the result given in equation (11) implies a constraint on the
isotropic dimension 5 SME parameter of

$k^{(5)}_{(V)00}\leq 4.2\times 10^{-34}\ {\rm GeV}^{-1}.$ (12)

Finally, it should be noted that with the redshift dependence obtained from
Equations (7) and (8), any reasonable redshift for a GRB similar to GRB041219a
and showing detectable polarization will give a constraint on $|\xi|$ below
$\sim 5\times 10^{-15}$ corresponding to a constraint on $k^{(5)}_{(V)00}$
below $\sim 10^{-33}$ GeV-1. This can be seen from Figure 1.

Much better tests of birefringence can be performed by polarization
measurements at higher $\gamma$-ray energies. The technology for measuring
polarization in the 5 to 100 MeV energy range using gas filled detectors is
now being developed and tested [25]. Studies of cosmological sources such as a
GRBs at such energies can probe values of $|\xi|$ several orders of magnitude
smaller than is presently possible.

##  5 Frame Independent Constraint

The vector $n$ in the EFT model given by equation (1) leads to strictly
isotropic physics only in one special frame, usually taken to be the frame in
which the cosmic microwave background is isotropic. In other frames the
dispersion relation will have anisotropic components. This can be taken into
account by using the general SME formalism [17]. There are then 16 independent
$k^{(5)}_{(jm)}$ parameters that are weighted by spherical harmonic
coefficients according to their spin weight with respect to the line of sight
unit vector ${\bf n}$. For GRB041219a this leads to the frame-independent
constraint

$|\sum_{jm}{Y_{jm}}(37^{\circ},0^{\circ})k^{(5)}_{(V)jm}|\leq 1.2\times
10^{-34}\ {\rm GeV^{-1}}.$ (13)

##  6 Other constraints and Implications

The Lorentz violating dispersion relation (2) implies that the group velocity
of photons, $v_{g}=1\pm\xi p/M_{Pl}$, is energy dependent. This leads to an
energy dependent dispersion in the arrival time at Earth for photons spread
over a finite energy range originating in a distant source. The result
obtained from observations of the $\gamma$-ray energy-time profile by the
Fermi satellite for the burst GRB090510 gives a limit of $\xi<0.82$ [26].
Thus, the time of flight constraint from Fermi, while still significant
because it gives $\xi<1$, remains many orders of magnitude weaker than the
birefringence constraint. However, the Fermi constraint is independent of the
EFT assumption of helicity dependence of the group velocity. Perhaps the best
constraint on LIV in general comes from a study of the highest energy cosmic
rays, giving a limit of $4.5\times 10^{-23}$ in the hadronic sector [5].

Thus, all of the present astrophysical data point to the conclusion that LIV
does not occur at the level $\xi(E/M_{Pl})$ with $\xi$ = $\cal{O}$(1). In
fact, in appears that $\xi\ll 1$. What this is telling us about the natures of
space-time and gravity at the Planck scale is still an open question.

## Acknowledgements

I would like to thank Neil Gehrels, Stanley Hunter, Sean Scully, Takanori
Sakamoto, Tonia Venters, and an anonymous referee for helpful discussions.

## References

  * [1] S. R. Coleman and S. L. Glashow, Phys. Rev. D 59, 116008 (1999).
  * [2] F. .W. Stecker and S. L. Glashow, Astropart. Phys 16, 97 (2001).
  * [3] M. Planck, Mitt. Thermodynamik, Folg. 5 (1899)
  * [4] T. Jacobson, S. Liberati and D. Mattingly, Annals of Physics 321, 150 (2006).
  * [5] F. W. Stecker and S. T. Scully, New J. Phys. 11 085003 (2009).
  * [6] R. C. Myers and M. Pospelov, Phys. Rev. Lett. 90, 211601 (2003).
  * [7] D. Colladay and V. A. Kostelecky, Phys. Rev. D 58, 116002 (1998).
  * [8] J. Ellis, N. E. Mavromatos, D. V. Nanopoulos, arXiv0912.3428 (2009) and references therein.
  * [9] D. Larson et al., Astrophys. J. Suppl. 192, 16 (2011)
  * [10] R. J. Gleiser and C. N. Kozameh, Phys. Rev. D 64, 083007 (2001).
  * [11] Y.-Z. Fan, D.-M. Wei and X. Dong, Mon. Not. Roy. Astr. Soc. 376, 1857 (2007).
  * [12] T. Jacobson, S. Liberati, D. Mattingly and F. Stecker, Phys. Rev. Lett. 93, 021101 (2004).
  * [13] W. Coburn and S. E. Boggs, Nature 423, 415 (2003).
  * [14] http://hesperia.gsfc.nasa.gov/hessi/
  * [15] R. E. Rutledge and D. B. Fox, Mon. Not. Roy. Astr. Soc. 350, 1288 (2004).
  * [16] C. Wigger et al., Astrophys. J. 613, 1088 (2004).
  * [17] V. A. Kostelecký and M. Mewes, Phys. Rev. D 80, 015020 (2009).
  * [18] V. A. Kostelecký and N. Russell, Rev. Mod. Phys. 83, 11 (2011).
  * [19] L. Maccione et al., Phys. Rev. D 78, 103003 (2008).
  * [20] S. Mc Glynn et al., Astron. and Astrophys. 466 895 (2007).
  * [21] D. Band et al., Astrophys. J. 413, 281 (1993).
  * [22] F.-Y. Wang, S. Qi and Z.-G. Dai, Mon. Not. Roy. Astr. Soc., in press (2011), arXiv:1105.0046.
  * [23] L. Xiao and B. E. Schaefer, Astrophys. J. 707, 387 (2009).
  * [24] S. McBreen et al. Astron. and Astrophys. 455 433 (2006).
  * [25] S. D. Hunter et al., in Space Telescopes and Instrumentation 1020: Ultraviolet to Gamma Rays, Proc. SPIE 7732, ed. et al., 773221 (2010).
  * [26] A. Abdo et al., Nature 462, 331 (2009).


# The number of perfect matchings, and the sandwich conjectures, of random
regular graphs

Pu Gao
University of Waterloo
<EMAIL_ADDRESS>
###### Abstract

We prove that the number of perfect matchings in ${\mathcal{G}}(n,d)$ is
asymptotically normal when $n$ is even, $d\to\infty$ as $n\to\infty$, and
$d=O(n^{1/7}/\log^{2}n)$. This is the first distributional result of spanning
subgraphs of ${\mathcal{G}}(n,d)$ when $d\to\infty$.

Moreover, we prove that ${\mathcal{G}}(n,d-1)$ and ${\mathcal{G}}(n,d)$ can be
coupled so that ${\mathcal{G}}(n,d-1)$ is a subgraph of ${\mathcal{G}}(n,d)$
with high probability when $d\to\infty$ and $d=o(n^{1/3})$. Further, if
$d=\Omega(\log^{7}n)$, $d=O(n^{1/7}/\log^{2}n)$, and $d\leq d^{\prime}\leq
n-1$ then ${\mathcal{G}}(n,d)$ and ${\mathcal{G}}(n,d^{\prime})$ can be
coupled so that asymptotically almost surely (a.a.s.) ${\mathcal{G}}(n,d)$ is
a subgraph of ${\mathcal{G}}(n,d^{\prime})$.

##  1 Introduction

We address two problems in this paper regarding ${\mathcal{G}}(n,d)$, the
random $d$-regular graph: the limiting distribution of the number of perfect
matchings in ${\mathcal{G}}(n,d)$, and the sandwich conjectures of
${\mathcal{G}}(n,d)$.

###  1.1 The number of perfect matchings

The study of subgraphs lies in the centre of random graph theory. The commonly
studied examples include spanning subgraphs such as perfect matchings,
Hamilton cycles, spanning trees, $H$-factors where $H$ has a fixed size, as
well as smaller subgraphs such as independent sets, cycles, and in general
subgraphs isomorphic to some given $H$ of fixed size. Let $Z_{H}$ denote the
number of subgraphs isomorphic to $H$. The phase transition of positive
$Z_{H}$ and the distribution of $Z_{H}$ have been well studied in
${\mathcal{G}}(n,p)$ and in ${\mathcal{G}}(n,m)$ for both small and large $H$.
It is interesting that $Z_{H}$ has different types of distributions for small
and large $H$. If $H$ has fixed size and is balanced, then $Z_{H}$ is
asymptotically normally distributed in ${\mathcal{G}}(n,p)$ and
${\mathcal{G}}(n,m)$ when ${\mathbb{E}}Z_{H}\to\infty$ [21]. The distribution
of $Z_{H}$ for $H$ with linear size becomes complicated. For $p\gg n^{-1/2}$
and $p$ not too close to 1, the numbers of perfect matchings, Hamilton cycles,
and spanning trees are log-normally distributed in ${\mathcal{G}}(n,p)$, but
are normally distributed in ${\mathcal{G}}(n,m)$ when $m\gg n^{3/2}$ [11]. For
$m=\Theta(n^{3/2})$, these random variables also become log-normally
distributed in ${\mathcal{G}}(n,m)$. It is not known if they remain log-
normally distributed in ${\mathcal{G}}(n,p)$ when $p=O(n^{-1/2})$ and in
${\mathcal{G}}(n,m)$ when $m=o(n^{3/2})$, although it is conjectured so [11].
The distributional phase transition of $Z_{H}$ when the size of $H$ grows from
constant to linear size has been studied in [7], where $H$ is the number of
$\ell$-matchings (i.e. matchings of size $\ell$). Its distribution in
${\mathcal{G}}(n,p)$ changes from normal to log-normal at the critical value
$\ell=\ell(p)\approx n\sqrt{p}$. Such distributional phase transition is also
observed in ${\mathcal{G}}(n,d)$ when $d$ is a fixed constant. It is well
known that the distributions of short cycles in ${\mathcal{G}}(n,d)$ are
asymptotically Poisson [1, 22], whereas the distribution of the number of
large subgraphs such as perfect matchings and Hamilton cycles in
${\mathcal{G}}(n,d)$ is of an unusual type [13] as follows. Suppose that $Z$
is the number of perfect matchings (or the number of Hamilton cycles) in
${\mathcal{G}}(n,d)$. Then the limiting distribution of the logarithm of
$Z/{\mathbb{E}}Z$ is an infinite linear combination of independent Poisson
variables. More precisely,

$\frac{Z}{{\mathbb{E}}Z}\to\prod_{i=1}^{\infty}(1+\delta_{i})^{X_{i}}e^{-\lambda_{i}\delta_{i}},\quad\mbox{as
$n\to\infty$},$ (1)

where $X_{1},X_{2},\ldots$ are independent Poisson variables with mean
$\lambda_{1},\lambda_{2},\ldots$, and $\delta_{1},\delta_{2},\ldots$ are real
numbers whose values depend on which subgraphs (i.e. perfect matchings or
Hamilton cycles) $Z$ counts. The distribution of $Z$ is determined by using
the small subgraph conditioning method, originally developed by Robinson and
Wormald [19, 20] to prove Hamiltonicity of ${\mathcal{G}}(n,d)$. The argument
is then tuned to produce the distribution result of $Z$ by Janson [13].
Recently, Greenhill, Isaev and Liang [10] proved that the number of spanning
trees in ${\mathcal{G}}(n,d)$ has the same type of distribution as (1). On the
other hand, Garmo [9] studied the distributional phase transition of the
number of $\ell$-cycles in ${\mathcal{G}}(n,d)$ as $\ell$ grows from constant
to linear in $n$. Its limiting distribution changes from a linear combination
of independent Poisson variables to the exponential of that form, and the
critical phase transition occurs when $\ell$ becomes linear in $n$.

It is natural to ask, in the case $d\to\infty$, whether the distribution type
of these subgraphs (e.g. perfect matchings, Hamilton cycles, spanning trees)
are the same as, or analogous to, that for constant $d$, and whether the
distributional phase transitions occur when the size of the subgraphs (e.g.
$\ell$-matchings and $\ell$-cycles) grows from constant to linear in $n$, as
for constant $d$. We give a negative answer to this question. There have been
few distributional results that are known for the number of subgraphs of
${\mathcal{G}}(n,d)$ when $d\to\infty$, even for small subgraphs. The limiting
distribution of the number of $\ell$-cycles was extended from constant $d$ and
$\ell$ to those such that $(d-1)^{2\ell-1}=o(n)$ by McKay, Wormald and Wysocka
[18]. Z. Gao and Wormald determined the limiting distributions of strictly
balanced graphs of fixed sizes for $d$ that grows sufficiently slowly with $n$
[8]. There has been no result on the distribution of the number of subgraphs
whose size is beyond $\log n$ when $d\to\infty$. In particular, the analysis
for the number of perfect matchings, Hamilton cycles, and the spanning trees
when $d=O(1)$, based on the configuration model, can no longer apply to obtain
the distributional result when $d\to\infty$.

One may expect that the number of large subgraphs such as perfect matchings or
Hamilton cycles would be of log-normal type in ${\mathcal{G}}(n,d)$ as
$d\to\infty$, which can be viewed as an analog of (1). One may also expect
that the number of $\ell$-matchings may exhibit a distributional phase
transition as $\ell$ grows from constant to linear in $n$, as what happens for
the $\ell$-matchings in ${\mathcal{G}}(n,p)$ and for the $\ell$-cycles in
${\mathcal{G}}(n,d)$ for constant $d$. In contrast with the intuition, we show
in this paper that the number of perfect matchings is asymptotically normally
distributed in ${\mathcal{G}}(n,d)$ when $d\to\infty$ as $n\to\infty$ and
$d=O(n^{1/7}/\log^{2}n)$. The power of the logarithmic term is not optimised.

######  Theorem 1.

Let $Y$ denote the number of perfect matchings in ${\mathcal{G}}(n,d)$ where
$n$ is even. Then $Y$ is asymptotically normally distributed if $d\to\infty$
as $n\to\infty$ and $d=O(n^{1/7}/\log^{2}n)$. More formally,

$\frac{Y-{\mathbb{E}}Y}{{\mathbb{E}}Y/\sqrt{6d^{3}}}\xrightarrow{d}{\mathcal{N}}(0,1),\quad\mbox{as
$n\to\infty$}.$

This result suggests that there is likely no distributional phase transition
on the number of $\ell$-matchings as $\ell$ grows. The condition
$d=O(n^{1/7}/\log^{2}n)$ in the result is imposed only for technical reasons
and we believe that the same distribution holds for all $d\to\infty$ until $d$
is too close to $n-1$.

To our knowledge this is the first result on the limiting distribution of the
number of spanning subgraphs in ${\mathcal{G}}(n,d)$ when $d\to\infty$. The
main contribution of Theroem 1 is the discovery of the distribution type of
the number of perfect matchings, and we believe that this phenomenon is
ubiquitous among other spanning subgraphs such as the number of Hamilton
cycles. For future research, it is interesting to determine the limiting
distributions of the number of $\ell$-matchings and $\ell$-cycles in
${\mathcal{G}}(n,d)$ for all $\ell$.

######  Conjecture 2.

The numbers of perfect matchings, Hamilton cycles, spanning trees, and
$k$-factors where $k\leq d-1$, are all asymptotically normally distributed in
${\mathcal{G}}(n,d)$ for all $d$ where $dn$ is even and
$\min\\{d,n-d\\}\to\infty$ as $n\to\infty$ (where $n$ is even in the case of
perfect matchings, and $kn$ is even in the case of $k$-factors).

######  Conjecture 3.

Suppose $\min\\{d,n-d\\}\to\infty$ as $n\to\infty$. The number of
$\ell$-cycles in ${\mathcal{G}}(n,d)$ is asymptotically normal for all
$3\leq\ell\leq n$. The number of $\ell$-matchings of ${\mathcal{G}}(n,d)$ is
asymptotically normal for all $3\leq\ell\leq n/2$.

######  Remark 4.

The condition $n-d\to\infty$ in the above conjectures are likely not
necessary. Indeed, when $n-d=o(n^{1/3})$ and $k=o(n^{1/3})$ the asymptotic
number of $k$-factors of a $d$-regular graph $G$ is independent of $G$ and can
be obtained by Theorem 13 in Section 3.

###  1.2 The sandwich conjectures of ${\mathcal{G}}(n,d)$

Analysis in ${\mathcal{G}}(n,d)$ is highly nontrivial, especially when
$d\to\infty$. Kim and Vu initiated the study of approximating
${\mathcal{G}}(n,d)$ by ${\mathcal{G}}(n,p)$, known as the sandwich conjecture
[14]. Since then the conjecture has been worked on by several groups [2, 5,
15, 3], and the conjecture is close to being fully resolved. Along the line of
the research there has been new conjectures that are proposed, one of which is
stated as follows [5, Conjecture 1.2].

######  Conjecture 5.

Let $0\leq d_{1}\leq d_{2}\leq n-1$ be integers, other than
$(d_{1},d_{2})=(1,2)$ or $(d_{1},d_{2})=(n-3,n-2)$. Assume that $d_{1}n$ and
$d_{2}n$ are both even. Then there exists a coupling $(G_{1},G_{2})$ such that
$G_{1}\sim{\mathcal{G}}(n,d_{1})$, $G_{2}\sim{\mathcal{G}}(n,d_{2})$ and
${\mathbb{P}}(G_{1}\subseteq G_{2})=1-o(1)$.

The conjecture is only known to be true for $(d_{1},d_{2})$ where $d_{1}=1$
and $3\leq d_{2}\leq n-1$, as well as for $(d_{1},d_{2})$ where $d_{2}-d_{1}$
is larger than some function of $d_{1}$ (see [5, Corollary 1.7] for the
precise statement). When $d_{1}$ and $d_{2}$ are both fixed constants and
$(d_{1},d_{2})\neq(1,2)$, it is known that ${\mathcal{G}}(n,d_{2})$ is
contiguous to the union of two independent copies of ${\mathcal{G}}(n,d_{1})$
and ${\mathcal{G}}(n,d_{2}-d_{1})$ conditional on ${\mathcal{G}}(n,d_{1})$ and
${\mathcal{G}}(n,d_{2}-d_{1})$ being disjoint (see more contiguity results in
[23, Section 4]). However, contiguity does not imply a coupling as in the
conjecture. In this paper we prove Conjecture 5 for a certain range of
$d_{1}$.

######  Theorem 6.

Conjecture 5 holds for all integers $d_{1}\leq d_{2}\leq n-1$ where
$d_{1}=\Omega(\log^{7}n)$ and $d_{1}=O(n^{1/7}/\log^{2}n)$ if $n$ is even.

Theorem 6 follows as a corollary of [3, Theorem 2] and the following theorem
that simultaneously couple a sequence of random regular graphs.

######  Theorem 7.

Suppose $d\to\infty$ and $d=O(n^{1/7}/\log^{2}n)$. For any
$\epsilon_{n}=o(1)$, there is a multiple coupling
$(G_{d},G_{d+1},\ldots,G_{\lfloor(1+\epsilon_{n})d\rfloor})$ such that
marginally $G_{i}\sim{\mathcal{G}}(n,i)$ for all $d\leq
i\leq\lfloor(1+\epsilon_{n})d\rfloor$ and jointly with probability $1-o(1)$,
$G_{d}\subseteq G_{d+1}\subseteq\cdots\subseteq
G_{\lfloor(1+\epsilon_{n})d\rfloor}$ simultaneously.

If we only consider $(d_{1},d_{2})$ where $d_{2}=d_{1}+1$ then we have the
following coupling theorem which holds for a much larger range of $d$.

######  Theorem 8.

Suppose $d\to\infty$ and $d=o(n^{1/3})$. There is a coupling $(G_{d},G_{d+1})$
where marginally $G_{d}\sim{\mathcal{G}}(n,d)$ and
$G_{d+1}\sim{\mathcal{G}}(n,d+1)$, and jointly with probability $1-o(1)$,
$G_{d}\subseteq G_{d+1}$.

Theorem 8 follows as a corollary of a more general version (Theorem 18) which
we state in Section 4. Theorem 8 indeed holds for all $d\to\infty$ and
$d=o(n^{1/2})$. However for a simpler proof we did not pursue that. See Remark
15 for more explanations.

##  2 Proof of Theroem 1

Recall that $Y$ denotes the number of perfect matchings in
${\mathcal{G}}(n,d)$. Throughout the paper we assume that $n$ is even. Let $X$
denote the number of triangles in ${\mathcal{G}}(n,d)$. We will approximate
$Y$ by a linear function of $X$ using linear regression, and then study the
distribution of $Y$ by the distribution of $X$. This method is known as
orthogonal decomposition and projection, developed by Janson [12]. Originally
this method is developed to determine the limiting distribution of the number
of (large) subgraphs in ${\mathcal{G}}(n,p)$, and Janson also applied the
method to determine the distributions of the numbers of spanning trees,
perfect matchings and Hamilton cycles in ${\mathcal{G}}(n,m)$ [11]. We are not
aware of its applications in ${\mathcal{G}}(n,d)$ before. More specifically,
we approximate $Y$ by $Y^{\star}=aX+b$ where $a={\bf Cov}(X,Y)/{\bf Var}X$ and
$b={\mathbb{E}}Y-a{\mathbb{E}}X$. The values of $a$ and $b$ are chosen so that
${\mathbb{E}}(Y-Y^{\star})=0$ and ${\mathbb{E}}((Y-Y^{\star})^{2})$ is
minimised. We prove that ${\mathbb{E}}((Y-Y^{\star})^{2})$ is sufficiently
small and thus the distribution of $Y$ is determined by the distribution of
$aX+b$. Since $X$ is asymptotically normally distributed, so is $Y$. The
expectation ${\mathbb{E}}X$, the variance ${\bf Var}X$ and the limiting
distribution of $X$ have been studied in [4, Theorems 8 and 10], which we
state below.

######  Theorem 9.

Suppose $d=o(n^{2/5})$. Then,

$\displaystyle{\mathbb{E}}X=\frac{(d-1)^{3}}{6}\left(1+O(1/n)\right),\quad{\bf
Var}X\sim{\mathbb{E}}X,$
$\displaystyle\frac{X-{\mathbb{E}}X}{\sqrt{{\mathbb{E}}X}}\xrightarrow{d}{\mathcal{N}}(0,1),\
\mbox{as $n\to\infty$}.$

Next, we calculate the expectation ${\mathbb{E}}Y$, the second moment
${\mathbb{E}}Y^{2}$ and the covariance ${\bf Cov}(X,Y)$, which allow us to
estimate $a$ and $b$ and to bound ${\mathbb{E}}((Y-Y^{\star})^{2})$.

######  Theorem 10.

$\displaystyle{\mathbb{E}}Y$
$\displaystyle=\frac{n!}{(n/2)!2^{n/2}}\left(\frac{e}{n}\right)^{n/2}\left(\frac{d-1}{d}\right)^{(\frac{d-1}{2})n}d^{\frac{n}{2}}\exp\left(\frac{1}{4}+O\left(\frac{d^{3}}{n}\right)\right)$
(2) $\displaystyle{\mathbb{E}}Y^{2}$
$\displaystyle=\left(1+\frac{1}{6d^{3}}+O\left(d^{-4}+\frac{d^{3}}{n}+\sqrt{\frac{d}{n}}\log^{6}n\right)\right)({\mathbb{E}}Y)^{2}$
(3) $\displaystyle{\bf Cov}(X,Y)$
$\displaystyle=\left(-\frac{1}{d^{3}}+O\left(d^{-4}+\frac{d}{n}\right)\right){\mathbb{E}}X{\mathbb{E}}Y.$
(4)

Proof of Theorem 1.  Recall the definition of $a$ and we make the following
claim.

######  Claim 11.

${\mathbb{E}}((Y-Y^{\star})^{2})=o(a^{2}\,{\bf Var}X)$.

Since

$\displaystyle Y$ $\displaystyle=$ $\displaystyle
Y^{\star}+(Y-Y^{\star})=a(X-{\mathbb{E}}X)+a{\mathbb{E}}X+b+(Y-Y^{\star}),$

it follows then that

$\displaystyle\frac{X-{\mathbb{E}}X}{\sqrt{{\mathbb{E}}X}}$ $\displaystyle=$
$\displaystyle\frac{Y-a{\mathbb{E}}X-b+(Y-Y^{\star})}{a\sqrt{{\mathbb{E}}X}}.$

By Claim 11, Theorem 9 and Markov’s inequality, asymptotically almost surely
(a.a.s.)

$|Y-Y^{\star}|=o(a\sqrt{{\bf Var}X})=o(a\sqrt{{\mathbb{E}}X}).$ (5)

Thus,

$\frac{Y-a{\mathbb{E}}X-b}{a\sqrt{{\mathbb{E}}X}}=\frac{X-{\mathbb{E}}X}{\sqrt{{\mathbb{E}}X}}-\frac{Y-Y^{\star}}{a\sqrt{{\mathbb{E}}X}},$

where the left hand side above is

$\frac{Y-{\mathbb{E}}Y}{a\sqrt{{\mathbb{E}}X}}\sim\frac{Y-{\mathbb{E}}Y}{{\mathbb{E}}Y/\sqrt{6d^{3}}},$

and the right hand side is a random variable whose limiting distribution is
${\mathcal{N}}(0,1)$ by (5) and Theorem 9. Consequently,

$\frac{Y-{\mathbb{E}}Y}{{\mathbb{E}}Y/\sqrt{6d^{3}}}\xrightarrow{d}{\mathcal{N}}(0,1),\quad\mbox{as
$n\to\infty$}.~{}~{}\vrule height=8.0pt,width=4.0pt,depth=0.0pt$

Proof of Claim 11.  By the definitions of $a$ and $b$,

$\displaystyle{\mathbb{E}}((Y-Y^{\star})^{2})$ $\displaystyle=$
$\displaystyle{\mathbb{E}}Y^{2}-2{\mathbb{E}}Y(aX+b)+{\mathbb{E}}(aX+b)^{2}$
$\displaystyle=$
$\displaystyle{\mathbb{E}}Y^{2}-2a{\mathbb{E}}XY-2({\mathbb{E}}Y-a{\mathbb{E}}X){\mathbb{E}}Y+a^{2}{\mathbb{E}}X^{2}+2a({\mathbb{E}}Y-a{\mathbb{E}}X){\mathbb{E}}X+({\mathbb{E}}Y)^{2}$
$\displaystyle-2a{\mathbb{E}}X{\mathbb{E}}Y+a^{2}({\mathbb{E}}X)^{2}$
$\displaystyle=$
$\displaystyle({\mathbb{E}}Y^{2}-({\mathbb{E}}Y)^{2})-2a({\mathbb{E}}XY-{\mathbb{E}}X{\mathbb{E}}Y)+a^{2}({\mathbb{E}}X^{2}-({\mathbb{E}}X)^{2})$
$\displaystyle=$ $\displaystyle{\bf Var}Y-2a{\bf Cov}(X,Y)+a^{2}{\bf Var}(X)$
$\displaystyle=$ $\displaystyle{\bf Var}Y-\frac{{\bf Cov}(X,Y)^{2}}{{\bf
Var}X}.$

By the definition of $a$, $a^{2}\,{\bf Var}X={\bf Cov}(X,Y)^{2}/{\bf Var}X$.
Hence it is sufficient to verify that

${\bf Var}Y-\frac{{\bf Cov}(X,Y)^{2}}{{\bf Var}X}\ll\frac{{\bf
Cov}(X,Y)^{2}}{{\bf Var}X}.$

It is thus sufficient to prove that

${\bf Var}Y\sim\frac{{\bf Cov}(X,Y)^{2}}{{\bf Var}X}.$ (6)

By Theorem 10,

${\bf
Var}Y={\mathbb{E}}Y^{2}-({\mathbb{E}}Y)^{2}=\left(\frac{1}{6d^{3}}+O\left(d^{-4}+\frac{d^{3}}{n}+\sqrt{\frac{d}{n}}\log^{6}n\right)\right)({\mathbb{E}}Y)^{2},$

and by Theorems 10 and 9,

$\frac{{\bf Cov}(X,Y)^{2}}{{\bf
Var}X}\sim\frac{\left(-\frac{1}{d^{3}}+O(d^{-4}+d/n)\right)^{2}({\mathbb{E}}X{\mathbb{E}}Y)^{2}}{{\mathbb{E}}X}\sim\frac{1}{6d^{3}}({\mathbb{E}}Y)^{2}.$

Now (6) follows, since $d\to\infty$ and $d=O(n^{1/7}/\log^{2}n)$ and thus
${\bf Var}Y\sim({\mathbb{E}}Y)^{2}/6d^{3}$.

##  3 Proof of Theorem 10

We will use the tools from [6, Theorem 1] and [17, Theorem 4.6] to estimate
${\mathbb{E}}Y$, ${\mathbb{E}}Y^{2}$ and ${\bf Cov}(X,Y)$.

###  3.1 Edge and subgraph probabilities in ${\mathcal{G}}(n,d)$

Let $H$ be a graph on $[n]$ and let ${\bf d}^{H}=(d_{1}^{H},\ldots,d_{n}^{H})$
denote the degree sequence of $H$. Suppose that $d_{i}^{H}\leq d$ for every
$1\leq i\leq n$. Let $H^{+}$ denote the event that $H$ is a subgraph of
${\mathcal{G}}(n,d)$. Let $|H|$ denote the number of edges in $H$. The
following result is a special case of [6, Theorem 1] for the conditional edge
probability ${\mathbb{P}}(uv\in{\mathcal{G}}(n,d)\mid H^{+})$.

######  Theorem 12.

Suppose $d=o(n)$ and suppose that $H$ is a graph on $[n]$ such that
$d_{i}^{H}\leq d$ for every $1\leq i\leq n$ and $dn-2|H|=\Omega(dn)$.

$\displaystyle{\mathbb{P}}\Big{(}uv\in{\mathcal{G}}(n,d)\mid
H^{+}\Big{)}=\left(1+O\left(\frac{d}{n}\right)\right)\frac{(d-d_{u}^{H})(d-d_{v}^{H})}{M-2|H|}.$

We will apply the following enumeration result of McKay [17, Theorem 4.6] to
estimate the probability of a (large) subgraph of ${\mathcal{G}}(n,d)$.

######  Theorem 13.

Let ${\bf g}=(g_{1},\ldots,g_{n})$ be a sequence of non-negative integers. Let
$m=m({\bf g})=\|{\bf g}\|_{1}/2$. Let $X$ be a simple graph on $[n]$ with
degree sequence ${\bf x}$. Let $\Delta({\bf g})$ and $\Delta({\bf x})$ denote
the maximum components of ${\bf g}$ and ${\bf x}$ respectively. Suppose
$\Delta({\bf g})\geq 1$, $\hat{\Delta}=o(m)$ where $\hat{\Delta}({\bf
g})=\Delta({\bf g})^{2}+\Delta({\bf g})\Delta({\bf x})$. Define

$\lambda=\lambda({\bf g})=\frac{1}{4m({\bf
g})}\sum_{j=1}^{n}(g_{i})_{2},\quad\mu=\mu({\bf g},X)=\frac{1}{2m({\bf
g})}\sum_{ij\in X}g_{i}g_{j}.$

Let $N({\bf g},X)$ denote the number of simple graphs with degree sequence
${\bf g}$ and with no edge in common with $X$. Then,

$N({\bf
g},X)=\frac{(2m)!}{m!2^{m}\prod_{j=1}^{n}g_{i}!}\exp\Big{(}-\lambda({\bf
g})-\lambda({\bf g})^{2}-\mu({\bf g},X)+O(\hat{\Delta}({\bf g})^{2}/m({\bf
g}))\Big{)}.$

######  Corollary 14.

Let $0\leq k\leq n/2$ be an integer. Let $H$ be a graph on $[n]$ containing
$k$ isolated edges and a collection of disjoint cycles spanning the remaining
$n-2k$ vertices. Then, with $\alpha=2k/n$,

$\displaystyle{\mathbb{P}}(H\subseteq{\mathcal{G}}(n,d))$
$\displaystyle=\frac{((d-2)n+2k)!\frac{dn}{2}!2^{n-k}d^{n}(d-1)^{n-2k}}{(\frac{(d-2)n}{2}+k)!(dn)!}\exp\left(\phi(d,\alpha)+O\left(\frac{d^{3}}{n}\right)\right)$
$\displaystyle=\left(\frac{e}{n}\right)^{(1-\frac{\alpha}{2})n}\left(\frac{d-2+\alpha}{d}\right)^{(\frac{d-2+\alpha}{2})n}d^{\frac{\alpha}{2}n}\left(d-1\right)^{(1-\alpha)n}\exp\left(\phi(d,\alpha)+O\left(\frac{d^{3}}{n}\right)\right),$

where

$\phi(d,\alpha)=\frac{4(d-2)^{2}-(d^{2}-5)\alpha^{2}-(2d^{2}-14d+20)\alpha}{4(d-2+\alpha)^{2}}.$
(7)

######  Remark 15.

Theorem 12 can be deduced from earlier work than [6], e.g. by McKay [16]. We
cite [6, Theorem 1] as it is written in form of conditional edge
probabilities, which is what we need in this paper.

A stronger version of Theorem 12 is available in [4, Theorem 6] which
estimates the conditional edge probabilities up to a relative error
$d^{2}/n^{2}$ instead of $d/n$. Using that result, we can deduce Corollary 14
with a smaller error $O(d^{2}/n)$ than $O(d^{3}/n)$. This will result in an
improvement in the range of $d$ in several of theorems in the paper, e.g. in
Theorems 8, 10, and 18. However, applying [4, Theorem 6] involves more
intensive calculations, and for a simpler proof we deduce Corollary 14 from
Theroem 13 instead.

Proof of Corollary 14.  Let ${\bf d}^{H}$ denote the degree sequence of $H$
and let ${\bf g}={\bf d}-{\bf d}^{H}$ where ${\bf d}=(d,\ldots,d)$. Then,
${\bf g}$ has exactly $2k$ components of value $d-1$ and $n-2k$ components of
value $d-2$. Hence,

$\displaystyle 2m({\bf g})$ $\displaystyle=(d-2)n+2k$
$\displaystyle\lambda({\bf g})$
$\displaystyle=\frac{1}{2((d-2)n+2k)}\Big{(}(d-1)(d-2)\cdot
2k+(d-2)(d-3)(n-2k)\Big{)}$ $\displaystyle\mu({\bf g},H)$
$\displaystyle=\frac{1}{(d-2)n+2k}\Big{(}(d-1)^{2}\cdot
k+(d-2)^{2}(n-2k)\Big{)}.$

and

$2m({\bf d})=dn,\quad\lambda({\bf
d})=\frac{1}{2dn}\left(d(d-1)n\right),\quad\mu({\bf d},\emptyset)=0.$

Moreover,

$\hat{\Delta}({\bf g}),\hat{\Delta}({\bf d})=O(d^{2})\quad\mbox{and}\quad
m({\bf g}),m({\bf d})=\Omega(dn).$

Thus,

$\displaystyle{\mathbb{P}}(H\subseteq{\mathcal{G}}(n,d))$
$\displaystyle=\frac{N({\bf g},H)}{N({\bf d},\emptyset)}$
$\displaystyle=\frac{((d-2)n+2k)!/\left(\frac{(d-2)n+2k}{2}!2^{\frac{(d-2)n+2k}{2}}\right)}{(dn)!/\left(\frac{dn}{2}!2^{\frac{dn}{2}}\right)}\cdot
d^{n}(d-1)^{n-2k}\exp(\phi(d,\alpha)+O(d^{3}/n)),$ (8)

where

$\displaystyle\phi(d,\alpha)$ $\displaystyle=-\lambda({\bf g})-\lambda({\bf
g})^{2}-\mu({\bf g},H)+\lambda({\bf d})+\lambda({\bf d})^{2}+\mu({\bf
d},\emptyset)$
$\displaystyle=\frac{4(d-2)^{2}-(d^{2}-5)\alpha^{2}-(2d^{2}-14d+20)\alpha}{4(d-2+\alpha)^{2}}.$
(9)

See Maple calculations of (9) in the Appendix. Now the corollary follows by
applying the Stirling formula to the factorials in (8). The relative error
$O(1/dn)$ in the Stirling formula is absorbed by $O(d^{3}/n)$.

###  3.2 ${\bf Cov}(X,Y)$

Fix a perfect matching $H$ of $K_{n}$. There are $n/2\cdot(n-2)$ ways to
choose a triangle $T$ such that $|H\cap T|=1$. By Theorem 12, the conditional
probability of $T$ given $H$ is $(d-1)^{2}(d-1)_{2}/(M-n)^{2}(1+O(d/n))$.
There are $\binom{n}{3}-(n/2\cdot(n-2))=(1+O(1/n))n^{3}/6$ ways to choose $T$
such that $H\cap T=\emptyset$. The conditional probability of $T$ given $H$ is
$(d-1)_{2}^{3}/(M-n)^{2}(1+O(d/n))$. Hence,

$\displaystyle{\mathbb{E}}XY$
$\displaystyle=\sum_{H}{\mathbb{P}}(H^{+})\sum_{T}{\mathbb{P}}(T\mid H^{+}),$

where the first summation is over all perfect matchings in $K_{n}$ and the
second summation is over all triangles in $K_{n}$. By the discussions above,
$\sum_{T}{\mathbb{P}}(T\mid H^{+})$ are the same for every perfect matching
$H$. Noting that $\sum_{H}{\mathbb{P}}(H^{+})={\mathbb{E}}Y$, we have

$\displaystyle{\mathbb{E}}XY$
$\displaystyle={\mathbb{E}}Y\left(\frac{n^{2}}{2}\frac{(d-1)^{3}(d-2)}{(M-n)^{2}}+\frac{n^{3}}{6}\cdot\frac{(d-1)^{3}(d-2)^{3}}{(M-n)^{3}}(1+O(d/n))\right)$
$\displaystyle={\mathbb{E}}Y\left(\frac{d^{2}}{2}(1-x)(1-2x)+\frac{d^{3}}{6}(1-2x)^{3}\right)(1+O(d/n))$
$\displaystyle=(1-1/d^{3}+O(1/d^{4}+d/n)){\mathbb{E}}X{\mathbb{E}}Y,$ (10)

where $x=1/d$ in the second equation above, and the last equation above is
obtained by taking the product of
$\frac{d^{2}}{2}(1-x)(1-2x)+\frac{d^{3}}{6}(1-2x)^{3}$ and
$({\mathbb{E}}X)^{-1}=(1+O(1/n))\frac{6}{d^{3}(1-x)^{3}}$ from Theorem 9, and
then taking the Taylor expansion of the product at $x=0$. We include the Maple
expansion formulae in the Appendix.

###  3.3 ${\mathbb{E}}Y$

Let $H$ be a perfect matching of $K_{n}$. By Corollary 14 with $k=n/2$ (i.e.
$\alpha=1$), we have $\phi(d,\alpha)=1/4$ and thus,

${\mathbb{P}}(H\subseteq{\mathcal{G}}(n,d))=(1+O(d^{3}/n))\rho_{1}(n,d).$

where

$\rho_{1}(n,d)=\left(\frac{e}{n}\right)^{n/2}\left(\frac{d-1}{d}\right)^{(\frac{d-1}{2})n}d^{\frac{n}{2}}\exp\left(\frac{1}{4}\right)$
(11)

Hence,

${\mathbb{E}}Y=\frac{n!}{(n/2)!2^{n/2}}\rho_{1}(n,d).$ (12)

###  3.4 ${\mathbb{E}}Y^{2}$

Let $0\leq k\leq n/2$ be an integer. Fix two perfect matchings $H_{1}$ and
$H_{2}$ of $K_{n}$ such that $|H_{1}\cap H_{2}|=k$. Let
$\alpha=\alpha(k)=2k/n$. Then, by Corollary 14,

${\mathbb{P}}(H_{1}\cup
H_{2}\subseteq{\mathcal{G}}(n,d))=(1+O(d^{3}/n))\rho_{2}(n,d,\alpha)$

where

$\displaystyle\rho_{2}(n,d,\alpha)$
$\displaystyle=\frac{((d-2)n+2k)!\frac{dn}{2}!2^{n-k}d^{n}(d-1)^{n-2k}}{(\frac{(d-2)n}{2}+k)!(dn)!}\exp\left(\phi(d,\alpha)\right)$
(13)
$\displaystyle=\left(\frac{e}{n}\right)^{(1-\frac{\alpha}{2})n}\left(\frac{d-2+\alpha}{d}\right)^{(\frac{d-2+\alpha}{2})n}d^{\frac{\alpha}{2}n}\left(d-1\right)^{(1-\alpha)n}\exp\left(\phi(d,\alpha)+O(1/dn)\right),$
(14)

with $\phi(d,\alpha)$ defined in (7). Next we count pairs $(H_{1},H_{2})$ of
perfect matchings of $K_{n}$ such that $|H_{1}\cap H_{2}|=k$.

######  Lemma 16.

The number of pairs $(H_{1},H_{2})$ of perfect matchings of $K_{n}$ such that
$|H_{1}\cap H_{2}|=k$ is

$(1+O((n-2k)^{-1}))\frac{n!}{2^{k}k!\sqrt{e\pi(n-2k)/2}}.$

Proof.  The exponential generating function for alternating cycles of length
at least 4 is

$F(z)=\sum_{n=2}^{\infty}\frac{(2n)!}{2\cdot 2n}\cdot
2\cdot\frac{z^{2n}}{(2n)!}=-\frac{1}{2}(\log(1-x^{2})+x^{2}).$

Thus, the number of pairs $(H_{1},H_{2})$ of disjoint perfect matchings of
$K_{2m}$ is

$(2m)!\cdot[z^{2m}]e^{F(z)}=(2m)!\cdot[z^{2m}]\frac{e^{-z^{2}/2}}{\sqrt{1-z^{2}}}=\frac{(2m)!}{\sqrt{e\pi
m}}(1+O(m^{-1})).$

Thus, the number of pairs $(H_{1},H_{2})$ of perfect matchings of $K_{n}$ such
that $|H_{1}\cap H_{2}|=k$ is

$\binom{n}{2k}\cdot\frac{(2k)!}{2^{k}k!}\cdot\frac{(n-2k)!}{\sqrt{e\pi(n-2k)/2}}(1+O((n-2k)^{-1}))=(1+O((n-2k)^{-1}))\frac{n!}{2^{k}k!\sqrt{e\pi(n-2k)/2}}.~{}~{}\vrule
height=8.0pt,width=4.0pt,depth=0.0pt$

By Lemma 16

$\displaystyle{\mathbb{E}}Y^{2}$
$\displaystyle=\sum_{k=0}^{n/2-1}(1+O((n-2k)^{-1}+d^{2}/n))\frac{n!}{2^{k}k!\sqrt{e\pi(n-2k)/2}}\rho_{2}(n,d,\alpha(k))+{\mathbb{E}}Y.$
(15)

Next, we show that the main contribution to ${\mathbb{E}}Y^{2}$ are from $k$
near some specific value. The proof of the lemma is postponed till Section
3.4.2.

######  Lemma 17.

Let

$\bar{\alpha}=\frac{1}{d},\quad\bar{k}=\lfloor\bar{\alpha}n/2\rfloor,\quad\bar{\delta}=\frac{2d}{n}\frac{d(d-2)}{(d-1)^{2}}.$

Then,

${\mathbb{E}}Y^{2}=(1+O(d^{3}/n))\frac{2}{\sqrt{e\bar{\delta}}}\frac{n!\rho_{2}(n,d,\bar{\alpha})}{\bar{k}!2^{\bar{k}}\sqrt{n-2\bar{k}}}.$

####  3.4.1 Comparing ${\mathbb{E}}Y^{2}$ with $({\mathbb{E}}Y)^{2}$

We complete the proof of Theorem 10 by verifying that

${\mathbb{E}}Y^{2}=\left(1+\frac{1}{6d^{3}}+O(\xi)\right)({\mathbb{E}}Y)^{2},$
(16)

where $\xi=d^{-4}+\frac{d^{3}}{n}+\sqrt{\frac{d}{n}}\log^{6}n$. By (12) and
Lemma 17,

$({\mathbb{E}}Y)^{2}=(1+O(\xi))2\left(\frac{n}{e}\right)^{n}\rho_{1}(n,d)^{2},$

and

${\mathbb{E}}Y^{2}=(1+O(\xi))\frac{2}{\sqrt{e\bar{\delta}}}\sqrt{\frac{n}{\bar{k}}}\frac{(n/e)^{n}\rho_{2}(n,d,\bar{\alpha})}{(2\bar{k}/e)^{\bar{k}}\sqrt{n-2\bar{k}}}.$

Hence,

$\displaystyle\frac{{\mathbb{E}}Y^{2}}{({\mathbb{E}}Y)^{2}}=(1+O(\xi))\sqrt{\frac{n}{e\bar{\delta}\bar{k}(n-2\bar{k})}}\frac{\rho_{2}(n,d,\bar{\alpha})}{(2\bar{k}/e)^{\bar{k}}\rho_{1}(n,d)^{2}}$

With straightforward but tedious calculations (see Appendix for more details)

$\sqrt{\frac{n}{e\bar{\delta}\bar{k}(n-2\bar{k})}}\frac{\rho_{2}(n,d,\bar{\alpha})}{(2\bar{k}/e)^{\bar{k}}\rho_{1}(n,d)^{2}}=\left(1+\frac{1}{6d^{3}}+O(\xi)\right),$
(17)

and now (16) follows.

####  3.4.2 Proof of Lemma 17

Recall from (18) that

$\displaystyle{\mathbb{E}}Y^{2}$
$\displaystyle=\sum_{k=0}^{n/2-1}(1+O((n-2k)^{-1}+d^{2}/n))\frac{n!}{2^{k}k!\sqrt{e\pi(n-2k)/2}}\rho_{2}(n,d,\alpha(k))+{\mathbb{E}}Y$
$\displaystyle=\sum_{k=0}^{n/2-1}(1+O((n-2k)^{-1}+d^{2}/n))\sqrt{\frac{2}{e\pi}}n!\varphi(k)+{\mathbb{E}}Y,$
(18)

where

$\varphi(k)=\frac{\rho_{2}(n,d,\alpha(k))}{2^{k}k!\sqrt{n-2k}}.$

The proof of Lemma 17 is standard. We prove that the summand in (18) is
maximised at $\bar{k}$. Then we approximate the summation around $\bar{k}$ by
an integral of a function of form $e^{-x^{2}}$. The contributions to (18) from
$k$ far away from $\bar{k}$ is negligible.

It is easy to see then that
$\exp(\phi(d,\alpha(k)-\phi(d,\alpha(k-1))))=\exp(O(1/n))$. Hence, by (13),

$\frac{\varphi(k)}{\varphi(k-1)}=\frac{(d-2)n+2k}{2(d-1)^{2}k}\left(1+O\left(\frac{1}{n}\right)\right)\quad\mbox{for
all $k$}.$ (19)

The maximum is achieved at $k=\bar{k}$. Moreover,

$\displaystyle\frac{\varphi(\bar{k}+j)}{\varphi(\bar{k})}$
$\displaystyle=\left(1+O\left(\frac{j}{n}\right)\right)\prod_{i=1}^{j}\frac{(d-2)n+2\bar{k}+2i}{2(\bar{k}+i)(d-1)^{2}}$
$\displaystyle=\left(1+O\left(\frac{j}{n}\right)\right)\prod_{i=1}^{j}\left(1+i\left(\frac{2}{(d-2)n+2\bar{k}}-\frac{1}{\bar{k}}\right)+O\left(\frac{i^{2}}{\bar{k}^{2}}\right)\right)$
$\displaystyle=\exp\left(-\bar{\delta}\sum_{i=1}^{j}i+O\left(\frac{j^{3}}{\bar{k}^{2}}+\frac{j}{n}\right)\right)$
$\displaystyle=\exp\left(-\frac{\bar{\delta}}{2}j^{2}+O\left(\bar{\delta}j+\frac{j^{3}}{\bar{k}^{2}}\right)\right),$
(20)

recalling that $\bar{\alpha}=1/d$, $\bar{k}=\lfloor\bar{\alpha}\cdot
n/2\rfloor$ and

$\displaystyle\bar{\delta}$
$\displaystyle=\frac{2d}{n}\frac{d(d-2)}{(d-1)^{2}}=\left(\frac{1}{\bar{k}}-\frac{2}{(d-2)n+2\bar{k}}\right)\left(1+O\left(\frac{1}{\bar{k}}\right)\right)$
$\displaystyle=\left(\frac{1}{\bar{k}}-\frac{2}{(d-2)n+2\bar{k}}\right)\left(1+O\left(\frac{d}{n}\right)\right).$

It follows then that

$\displaystyle\sum_{\bar{k}-\bar{\delta}^{-1/2}\log^{2}(1/\bar{\delta})\leq
k\leq\bar{k}+\bar{\delta}^{-1/2}\log^{2}(1/\bar{\delta})}\varphi(k)$
$\displaystyle=\left(1+O\left(\xi\right)\right)\varphi(\bar{k})\sqrt{\frac{2}{\bar{\delta}}}\int_{-\infty}^{+\infty}e^{-x^{2}}dx=\left(1+O\left(\xi\right)\right)\sqrt{\frac{2\pi}{\bar{\delta}}}\varphi(\bar{k}),$

where

$\xi=\bar{\delta}^{1/2}\log^{2}(1/\bar{\delta})+\frac{\bar{\delta}^{-3/2}\log^{6}(1/\bar{\delta})}{\bar{k}^{2}}+\bar{\delta}^{1/2}=O\left(\sqrt{\frac{d}{n}}\log^{6}n\right).$

Note that the first two terms in $\xi$ comes from the accumulative error
$O(\bar{\delta}j+j^{3}/\bar{k}^{2})$ in (20), and the last term comes from
approximating the sum of $\exp(-\bar{\delta}j^{2}/2)$ by an integral. The
contributions to (18) from $k$ where
$|k-\bar{k}|>\bar{\delta}^{-1/2}\log^{2}(1/\bar{\delta})$ is smaller than
$n^{-1}$ as a relative error by (19) and standard calculations on summing a
geometrically bounded series. So

${\mathbb{E}}Y^{2}=\left(1+O\left(\sqrt{\frac{d}{n}}\log^{6}n+\frac{d^{2}}{n}\right)\right)\frac{2}{\sqrt{e\bar{\delta}}}\frac{n!\rho_{2}(n,d,\bar{\alpha})}{\bar{k}!2^{\bar{k}}\sqrt{n-2\bar{k}}}.~{}~{}\vrule
height=8.0pt,width=4.0pt,depth=0.0pt$

##  4 Proofs of Theorems 6 – 8

###  4.1 Proof of Theorem 6

We prove Theorem 6 assuming Theorem 7. Let $d_{1}=\Omega(\log^{7}n)$. By [3,
Theorem 2], there exists $\delta_{n}=o(1)$ such that ${\mathcal{G}}(n,d_{1})$
can be coupled with ${\mathcal{G}}(n,(1+\delta_{n})d_{1}/n)$ such that a.a.s.
${\mathcal{G}}(n,d_{1})\subseteq{\mathcal{G}}(n,(1+\delta_{n})d_{1}/n)$.
Suppose $d_{2}\geq(1+2\delta_{n})d_{1}$. Then, again by [3, Theorem 2], a.a.s.
${\mathcal{G}}(n,d_{2})$ can be coupled with
${\mathcal{G}}(n,(1+\delta_{n})d_{1}/n)$ such that a.a.s.
${\mathcal{G}}(n,(1+\delta_{n})d_{1}/n)\subseteq{\mathcal{G}}(n,d_{2})$. It
follows now that Conjecture 5 holds for any $(d_{1},d_{2})$ where
$d_{1}=\Omega(\log^{7}n)$ and $d_{2}-d_{1}\geq 2\delta_{n}d_{1}$. Now suppose
$d_{1}+1\leq d_{2}<(1+2\delta_{n})d_{1}$ and $d_{1}=O(n^{1/7}/\log^{2}n)$.
Then, there is a coupling where a.a.s.
${\mathcal{G}}(n,d_{1})\subseteq{\mathcal{G}}(n,d_{2})$ by Theorem 7. Now
Theroem 6 follows.

###  4.2 Couple ${\mathcal{G}}(n,d)$ and ${\mathcal{G}}(n,d+1)$

Throughout this section we assume $d\to\infty$ and $d=o(n^{1/3})$. Our goal is
to couple ${\mathcal{G}}(n,d)$ with ${\mathcal{G}}(n,d+1)$ so that
${\mathcal{G}}(n,d)\subseteq{\mathcal{G}}(n,d+1)$ with sufficiently high
probability. In the next section, we “stitch” a sequence of such couplings
together to obtain a simultaneous coupling as in Theorem 7.

Given $\alpha=\alpha_{n}=o(1)$, define $\eta=\eta(\alpha)$ where

$\eta(\alpha)=2\alpha+\frac{1}{d^{3}\alpha^{2}}+\frac{C^{\prime}d^{3}}{n\alpha^{2}}+\frac{C^{\prime}\sqrt{d/n}\log^{6}n}{\alpha^{2}},\quad\mbox{where
$C^{\prime}>0$ is a sufficiently large constant.}$

We prove the following stronger version of Theorem 8.

######  Theorem 18.

Assume $\alpha=o(1)$ is such that $\eta(\alpha)=o(1)$. There is a coupling
$(G_{d},G_{d+1})$ where marginally $G_{d}\sim{\mathcal{G}}(n,d)$ and
$G_{d+1}\sim{\mathcal{G}}(n,d+1)$, and jointly with probability $1-O(\eta)$,
$G_{d}\subseteq G_{d+1}$.

####  4.2.1 The coupling procedure

For $G\in{\mathcal{G}}(n,d+1)$ let $Y(G)$ denote the number of perfect
matchings of $G$. For $G\in{\mathcal{G}}(n,d)$ let $Z(G)$ denote the number of
perfect matchings in $K_{n}\setminus G$. We say $G$ and $G^{\prime}$ are
related, denoted by $G\sim G^{\prime}$, for $G\in{\mathcal{G}}(n,d)$ and
$G^{\prime}\in{\mathcal{G}}(n,d+1)$, if $G\subseteq G^{\prime}$. We can
represent this relation using an auxiliary directed bipartite graph $X$ where
$V(X)={\mathcal{G}}(n,d)\cup{\mathcal{G}}(n,d+1)$, and $(G,G^{\prime})$ is an
arc if $G\subseteq G^{\prime}$. Thus, a $d$-regular graph $G$ in $X$ has out-
degree $Z(G)$, and a $(d+1)$-regular graph $G^{\prime}$ in $X$ has in-degree
$Y(G^{\prime})$.

Let $\alpha=\alpha_{n}=o(1)$. By Theorem 10 and Chebyshev’s inequality,

${\mathbb{P}}_{{\mathcal{G}}(n,d+1)}\left(|Y-{\mathbb{E}}Y|\geq\alpha{\mathbb{E}}Y\right)\leq\frac{{\bf
Var}Y}{\alpha^{2}({\mathbb{E}}Y)^{2}}=\frac{1}{6d^{3}\alpha^{2}}+O\left(\frac{1}{d^{4}\alpha^{2}}+\frac{d^{3}}{n\alpha^{2}}+\frac{\sqrt{d/n}\log^{6}n}{\alpha^{2}}\right).$
(21)

By Theorem 13 (wtih ${\bf g}$ being the all one vector and $X$ being a
$d$-regular graph) and Theorem 10, there exists a sufficiently large constant
$C>0$ such that

$\displaystyle|Z(G)-Z^{\star}|$ $\displaystyle\leq C\frac{d^{2}}{n}\cdot
Z^{\star},\ \mbox{for all $d$-regular graph $G$},$ (22)
$\displaystyle|{\mathbb{E}}Y(G)-Y^{\star}|$ $\displaystyle\leq
C\frac{d^{3}}{n}\cdot Y^{\star},\ \mbox{for $G\sim{\mathcal{G}}(n,d+1)$},$
(23)

where

$Z^{\star}=\frac{n!e^{-d/2}}{(n/2)!2^{n/2}},\quad
Y^{\star}=\frac{n!e^{1/4}}{(n/2)!2^{n/2}}\left(\frac{e}{n}\right)^{n/2}\left(\frac{d-1}{d}\right)^{(\frac{d-1}{2})n}d^{\frac{n}{2}}$

Let

$\displaystyle\underline{Y}$ $\displaystyle=(1-\alpha-Cd^{3}/n)\,Y^{\star}$
(24) $\displaystyle\overline{Z}$
$\displaystyle=\left(1+C\frac{d^{2}}{n}\right)Z^{\star}$ (25)

Define

$\displaystyle{\mathcal{B}}$ $\displaystyle=\\{G\in{\mathcal{G}}(n,d+1):\
Y(G)<\underline{Y}\\}$ (26) $\displaystyle{\mathcal{B}}^{\prime}$
$\displaystyle=\\{G\in{\mathcal{G}}(n,d+1):\
Y(G)>(1+\alpha+Cd^{3}/n)\,Y^{\star}\\}.$ (27)

By (21) and (23),

${\mathbb{P}}({\mathcal{B}}\cup{\mathcal{B}}^{\prime})=\frac{1}{6d^{3}\alpha^{2}}+O\left(\frac{1}{d^{4}\alpha^{2}}+\frac{d^{3}}{n\alpha^{2}}+\frac{\sqrt{d/n}\log^{6}n}{\alpha^{2}}\right).$
(28)

Let $D$ and $\hat{D}$ denote the total in-degrees of ${\mathcal{G}}(n,d+1)$
and ${\mathcal{G}}(n,d+1)\setminus{\mathcal{B}}$ respectively in $X$. That is,

$\displaystyle D$
$\displaystyle=|\\{(G,G^{\prime})\in{\mathcal{G}}(n,d)\times{\mathcal{G}}(n,d+1):G\sim
G^{\prime}\\}|$ (29) $\displaystyle\hat{D}$
$\displaystyle=|\\{(G,G^{\prime})\in{\mathcal{G}}(n,d)\times({\mathcal{G}}(n,d+1)\setminus{\mathcal{B}}):G\sim
G^{\prime}\\}|.$ (30)

Further, let $d^{-}({\mathcal{B}})$ and $d^{-}({\mathcal{B}}^{\prime})$ denote
the total in-degrees of ${\mathcal{B}}$ and ${\mathcal{B}}^{\prime}$
respectively in $X$. We prove the following bounds on $d^{-}({\mathcal{B}})$
and $d^{-}({\mathcal{B}}^{\prime})$. Recall that

$\eta(\alpha)=2\alpha+\frac{1}{d^{3}\alpha^{2}}+\frac{C^{\prime}d^{3}}{n\alpha^{2}}+\frac{C^{\prime}\sqrt{d/n}\log^{6}n}{\alpha^{2}},$

where $C^{\prime}>0$ is a sufficiently large constant.

######  Lemma 19.

Assume $\alpha=o(1)$ is such that $\eta(\alpha)=o(1)$. Then,
$d^{-}({\mathcal{B}})+d^{-}({\mathcal{B}}^{\prime})\leq\eta D$.

Proof.  The number of edges in $X$ is $D=|{\mathcal{G}}(n,d+1)|{\mathbb{E}}Y.$
This can be rewritten as

$d^{-}({\mathcal{B}})+d^{-}({\mathcal{B}}^{\prime})+|{\mathcal{G}}(n,d+1)\setminus(B\cup
B^{\prime})|\cdot{\mathbb{E}}Y\left(1+\xi\right),$

where $|\xi|\leq\alpha+O(d^{3}/n)$ since
$|Y(G)/{\mathbb{E}}Y-1|\leq\alpha+O(d^{3}/n)$ for all
$G\notin{\mathcal{B}}\cup{\mathcal{B}}^{\prime}$ by the definition of
${\mathcal{B}}$ and ${\mathcal{B}}^{\prime}$. By (28), the above is equal to

$d^{-}({\mathcal{B}})+d^{-}({\mathcal{B}}^{\prime})+|{\mathcal{G}}(n,d+1)|\left(1-\xi^{\prime}\right)\cdot{\mathbb{E}}Y\left(1+\xi\right)$

where

$0\leq\xi^{\prime}\leq\frac{1}{6d^{3}\alpha^{2}}+O\left(\frac{1}{d^{4}\alpha^{2}}+\frac{d^{3}}{n\alpha^{2}}+\frac{\sqrt{d/n}\log^{6}n}{\alpha^{2}}\right).$

Thus,

$|{\mathcal{G}}(n,d+1)|{\mathbb{E}}Y=d^{-}({\mathcal{B}})+d^{-}({\mathcal{B}}^{\prime})+|{\mathcal{G}}(n,d+1)|\cdot{\mathbb{E}}Y(1-\xi^{\prime}+\xi-\xi^{\prime}\xi),$

which implies that

$d^{-}({\mathcal{B}})+d^{-}({\mathcal{B}}^{\prime})=|{\mathcal{G}}(n,d+1)|\cdot{\mathbb{E}}Y(\xi^{\prime}-\xi+\xi^{\prime}\xi)<\eta
D,$

by the definition of $\eta$.

Finally we are ready to define the coupling $(G_{d},G_{d+1})$.

  * •

Let $G_{d}$ be a uniformly random graph in ${\mathcal{G}}(n,d)$ and let
$G^{\prime}$ be the graph obtained from $G_{d}$ by adding a uniformly random
perfect matching of $K_{n}\setminus G_{d}$. Let $H$ be a uniformly random
graph in ${\mathcal{G}}(n,d+1)$ independent of $G^{\prime}$.

  * •

If $G^{\prime}\in{\mathcal{B}}$ then let $G_{d+1}=G^{\prime}$ with probability
$(1-\eta)\frac{Z(G)}{\overline{Z}}$ and let $G_{d+1}=H$ with the remaining
probability.

  * •

If $G^{\prime}\in{\mathcal{G}}(n,d+1)\setminus{\mathcal{B}}$, then

$G_{d+1}=\left\\{\begin{array}[]{ll}G^{\prime}&\mbox{with probability
$(1-\eta)\frac{Z(G)}{\overline{Z}}\cdot\frac{\underline{Y}}{Y(G^{\prime})}$}\\\
G^{\prime\prime}&\mbox{with probability
$(1-\eta)\frac{Z(G)}{\overline{Z}}\frac{(\underline{Y}-Y(G^{\prime\prime}))}{\hat{D}}$
for every $G^{\prime\prime}\in{\mathcal{B}}$}\\\ H&\mbox{with the remaining
probability}.\end{array}\right.$

The following lemma justifies that the coupling procedure is well defined.

######  Lemma 20.

$\displaystyle\frac{Z(G)}{\overline{Z}}$ $\displaystyle\geq
1-3C\frac{d^{2}}{n},\quad\mbox{for every $d$-regular graph $G$}$
$\displaystyle\frac{\underline{Y}}{Y(G^{\prime})}$ $\displaystyle\geq
1-3\alpha-3C\frac{d^{3}}{n},\quad\mbox{for every $(d+1)$-regular graph
$G\in{\mathcal{G}}(n,d+1)\setminus({\mathcal{B}}\cup{\mathcal{B}}^{\prime})$,}$

and for every
$(G,G^{\prime})\in{\mathcal{G}}(n,d)\times\big{(}{\mathcal{G}}(n,d+1)\setminus{\mathcal{B}}\big{)}$
where $G\sim G^{\prime}$,

$(1-\eta)\frac{Z(G)}{\overline{Z}}\cdot\frac{\underline{Y}}{Y(G^{\prime})}+\sum_{G^{\prime\prime}\in{\mathcal{B}}}(1-\eta)\frac{Z(G)}{\overline{Z}}\frac{(\underline{Y}-Y(G^{\prime\prime}))}{\hat{D}}\leq
1.$

Proof.  The first inequality in the lemma follows by (22) and (25). The second
inequality in the lemma follows by (24), (26) and (27). For the last
inequality, note that

$\frac{Z(G)}{\overline{Z}},\frac{\underline{Y}}{Y(G^{\prime})}\leq 1$

always and thus it is sufficient to show that

$\sum_{G^{\prime\prime}\in{\mathcal{B}}}\frac{(\underline{Y}-Y(G^{\prime\prime}))}{\hat{D}}\leq\eta.$

By (21),

$|{\mathcal{B}}|\leq\left(\frac{1}{6d^{3}\alpha^{2}}+O\left(\frac{1}{d^{4}\alpha^{2}}+\frac{d^{3}}{n\alpha^{2}}+\frac{\sqrt{d/n}\log^{6}n}{\alpha^{2}}\right)\right)|{\mathcal{G}}(n,d+1)|.$

By Lemma 19,

$\hat{D}=(1+O(\eta))D=(1+O(\eta))|{\mathcal{G}}(n,d+1)|{\mathbb{E}}Y=(1+O(\eta))|{\mathcal{G}}(n,d+1)|\underline{Y}.$

Thus,

$\sum_{G^{\prime\prime}\in{\mathcal{B}}}(\underline{Y}-Y(G^{\prime\prime}))\leq\underline{Y}|{\mathcal{B}}|\leq\left(\frac{1}{6d^{3}\alpha^{2}}+O\left(\frac{1}{d^{4}\alpha^{2}}+\frac{d^{3}}{n\alpha^{2}}+\frac{\sqrt{d/n}\log^{6}n}{\alpha^{2}}\right)\right)\cdot|{\mathcal{G}}(n,d+1)|\underline{Y}\leq\eta\hat{D},$

by the definition of $\eta$. Thus, the last inequality of the lemma follows.

####  4.2.2 Proof of Theorem 18

By the construction, $G_{d}$ is obviously distributed as ${\mathcal{G}}(n,d)$
marginally. We prove that the marginal distribution of $G_{d+1}$ is
${\mathcal{G}}(n,d+1)$. Let $G^{\prime}$ be a $(d+1)$-regular graph. Let

$\sigma_{d}=\frac{1}{|{\mathcal{G}}(n,d)|},\quad\sigma_{d+1}=\frac{1}{|{\mathcal{G}}(n,d+1)|}.$

If $\hat{G}\in{\mathcal{G}}(n,d+1)\setminus{\mathcal{B}}$ then

${\mathbb{P}}(G_{d+1}=\hat{G})=\sum_{G:G\sim\hat{G}}\frac{\sigma_{d}}{Z(G)}\cdot(1-\eta)\frac{Z(G)}{\overline{Z}}\cdot\frac{\underline{Y}}{Y(\hat{G})}+\varphi,$
(31)

where

$\displaystyle\varphi=$
$\displaystyle\sum_{\begin{subarray}{c}(G,G^{\prime}):\\\
G^{\prime}\in{\mathcal{B}}\end{subarray}}\frac{\sigma_{d}}{Z(G)}\left(1-(1-\eta)\frac{Z(G)}{\overline{Z}}\right)\sigma_{d+1}$
$\displaystyle+\sum_{\begin{subarray}{c}(G,G^{\prime}):\\\ G^{\prime}\notin
B\end{subarray}}\frac{\sigma_{d}}{Z(G)}\left(1-(1-\eta)\frac{Z(G)}{\overline{Z}}\cdot\frac{\underline{Y}}{Y(G^{\prime})}-\sum_{G^{\prime\prime}\in{\mathcal{B}}}(1-\eta)\frac{Z(G)}{\overline{Z}}\frac{(\underline{Y}-Y(G^{\prime\prime}))}{\hat{D}}\right)\sigma_{d+1}.$

In the first summation in (31), $\sigma_{d}/Z(G)$ is the probability that
$G_{d}=G$ and $G^{\prime}=\hat{G}$. Conditioning on that,
$(1-\eta)\frac{Z(G)}{\overline{Z}}\cdot\frac{\underline{Y}}{Y(G^{\prime})}$ is
the probability that $G_{d+1}=G^{\prime}$. Thus this summation gives the
contribution to ${\mathbb{P}}(G_{d+1}=\hat{G})$ from the case that
$G^{\prime}=\hat{G}$ and $G_{d+1}$ is set $G^{\prime}$. Similarly, it is easy
to see that $\varphi$ is the probability that $G_{d+1}$ is set as a uniformly
random graph in ${\mathcal{G}}(n,d+1)$, no matter $G^{\prime}$ is in
${\mathcal{B}}$ or not. Note that the value of $\varphi$ is independent of
$\hat{G}$. Hence, by noting that $Y(\hat{G})=|\\{G:G\sim\hat{G}\\}|$, we
obtain

${\mathbb{P}}(G_{d+1}=\hat{G})=(1-\eta)\sigma_{d}\frac{\underline{Y}}{\overline{Z}}+\varphi,$

which is independent of $\hat{G}$ for all
$\hat{G}\in{\mathcal{G}}(n,d+1)\setminus{\mathcal{B}}$.

Next, suppose $\hat{G}\in{\mathcal{B}}$. Then,

${\mathbb{P}}(G_{d+1}=\hat{G})=\sum_{G:G\sim\hat{G}}\frac{\sigma_{d}}{Z(G)}\cdot(1-\eta)\frac{Z(G)}{\overline{Z}}+\sum_{\begin{subarray}{c}(G,G^{\prime}):\\\
G^{\prime}\notin{\mathcal{B}}\end{subarray}}\frac{\sigma_{d}}{Z(G)}\cdot(1-\eta)\frac{Z(G)}{\overline{Z}}\frac{(\underline{Y}-Y(\hat{G}))}{\hat{D}}+\varphi,$

where the second summation above is from the case where $G_{d}=G$,
$G^{\prime}\notin{\mathcal{B}}$ and $G_{d+1}$ is set to be $\hat{G}$ which
occurs with probability
$(1-\eta)\frac{Z(G)}{\overline{Z}}\frac{(\underline{Y}-Y(\hat{G}))}{D}$, given
$(G,G^{\prime})$. The first summation above gives
$(1-\eta)\sigma_{d}Y(\hat{G})/\overline{Z}$. The second summation above gives
$(1-\eta)\sigma_{d}(\underline{Y}-Y(\hat{G}))/\overline{Z}$ by (30). Hence,

${\mathbb{P}}(G_{d+1}=\hat{G})=(1-\eta)\sigma_{d}\frac{\underline{Y}}{\overline{Z}}+\varphi,$

which is independent of $\hat{G}$ for all $\hat{G}\in{\mathcal{B}}$, and is
the same for all $\hat{G}\in{\mathcal{G}}(n,d+1)\setminus{\mathcal{B}}$. This
confirms that the marginal distribution of $G_{d+1}$ is uniform in
${\mathcal{G}}(n,d+1)$.

Finally, we prove that $G_{d}\subseteq G_{d+1}$ with probability $1-O(\eta)$.
Note that $G_{d}\subseteq G_{d+1}$ if $G_{d+1}=G^{\prime}$ in the
construction. Thus, it is sufficient to show that the probability that
$G_{d+1}=H$ or $G_{d+1}=G^{\prime\prime}$ for some
$G^{\prime\prime}\in{\mathcal{B}}$ in the case $G^{\prime}\notin{\mathcal{B}}$
is $O(\eta)$.

Suppose $G^{\prime}\in{\mathcal{B}}$.

$\displaystyle{\mathbb{P}}(G_{d+1}=H)$
$\displaystyle=\sum_{G\in{\mathcal{G}}(n,d)}\sum_{G^{\prime}:G^{\prime}\sim
G}\frac{\sigma_{d}}{Z(G)}\left(1-(1-\eta)\frac{Z(G)}{\overline{Z}}\right)$
$\displaystyle=\sum_{G\in{\mathcal{G}}(n,d)}\sum_{G^{\prime}:G^{\prime}\sim
G}\frac{\sigma_{d}}{Z(G)}\left(1-(1-\eta)\left(1-O\left(\frac{d^{2}}{n}\right)\right)\right)\quad\mbox{(by
Lemma~{}\ref{lem:bounds})}$
$\displaystyle=O\left(\eta\right)\sum_{G\in{\mathcal{G}}(n,d)}\sigma_{d}\sum_{G^{\prime}:G^{\prime}\sim
G}\frac{1}{Z(G)}=O(\eta).$

Suppose $G^{\prime}\notin{\mathcal{B}}$.

$\displaystyle{\mathbb{P}}\Big{(}G_{d+1}=H\ \mbox{or}\
G_{d+1}=G^{\prime\prime}\ \mbox{for some
$G^{\prime\prime}\in{\mathcal{B}}$}\Big{)}$ $\displaystyle\hskip
14.22636pt=\sum_{G\in{\mathcal{G}}(n,d)}\sum_{\begin{subarray}{c}G^{\prime}:G^{\prime}\sim
G\\\
G^{\prime}\notin{\mathcal{B}}\end{subarray}}\frac{\sigma_{d}}{Z(G)}\left(1-(1-\eta)\frac{Z(G)}{\overline{Z}}\cdot\frac{\underline{Y}}{Y(G^{\prime})}\right)$
$\displaystyle\hskip
14.22636pt\leq\sum_{G\in{\mathcal{G}}(n,d)}\sum_{\begin{subarray}{c}G^{\prime}:G^{\prime}\sim
G\\\
G^{\prime}\notin{\mathcal{B}}\cup{\mathcal{B}}^{\prime}\end{subarray}}\frac{\sigma_{d}}{Z(G)}O(\eta)+\sum_{G\in{\mathcal{G}}(n,d)}\sum_{\begin{subarray}{c}G^{\prime}:G^{\prime}\sim
G\\\
G^{\prime}\in{\mathcal{B}}^{\prime}\end{subarray}}\frac{\sigma_{d}}{Z(G)},$

as for every $G\notin{\mathcal{B}}\cup{\mathcal{B}}^{\prime}$,
$1-(1-\eta)\frac{Z(G)}{\overline{Z}}\cdot\frac{\underline{Y}}{Y(G^{\prime})}=O(\eta)$
by Lemma 20 and for $G\in{\mathcal{B}}\cup{\mathcal{B}}^{\prime}$ we use the
trivial upper bound
$1-(1-\eta)\frac{Z(G)}{\overline{Z}}\cdot\frac{\underline{Y}}{Y(G^{\prime})}\leq
1$. Since $|\\{G^{\prime}:G^{\prime}\sim
G,G^{\prime}\notin{\mathcal{B}}\cup{\mathcal{B}}^{\prime}\\}|\leq|\\{G^{\prime}:G^{\prime}\sim
G\\}|=Z(G)$, and $\sigma_{d}\cdot|{\mathcal{G}}(n,d)|=1$, the first double
summation above is $O(\eta)$. The second double summation above is equal to

$\displaystyle\left(1+O(d^{2}/n)\right)\frac{\sigma_{d}}{Z^{\star}}|\\{(G,G^{\prime})\in{\mathcal{G}}(n,d)\times{\mathcal{B}}^{\prime}:G\sim
G^{\prime}\\}|=\left(1+O(d^{2}/n)\right)\frac{\sigma_{d}}{Z^{\star}}d^{-}({\mathcal{B}}^{\prime})$
$\displaystyle\hskip
14.22636pt\leq\left(1+O(d^{2}/n)\right)\frac{\sigma_{d}}{Z^{\star}}\eta
D\quad\mbox{(by Lemma~{}\ref{lem:degree-B})}$ $\displaystyle\hskip
14.22636pt=\left(1+O(d^{2}/n)\right)\frac{\sigma_{d}}{Z^{\star}}\eta|{\mathcal{G}}(n,d)|{\mathbb{E}}Z=O(\eta),$

where the last equality above holds because $Z^{\star}\sim{\mathbb{E}}Z$ and
$\sigma_{d}|{\mathcal{G}}(n,d)|=1$.

###  4.3 Proof of Theorem 8

Suppose $d\to\infty$ and $d=o(n^{1/3})$. Then there exists $\alpha=o(1)$ such
that $\eta(\alpha)=o(1)$. Theorem 8 follows by Theorem 18 with such a choice
of $\alpha$.

###  4.4 Proof of Theorem 7.

Suppose $d\to\infty$ and $d=O(n^{1/7}/\log^{2}n)$. Let $\alpha=1/d$. It
follows now that $\eta=O(1/d)$. We prove that for each $1\leq
j\leq\lfloor\epsilon_{n}d\rfloor$, there is a constant $C>0$, and a coupling
$(G_{d},\ldots,G_{d+j})$ where $G_{i}\sim{\mathcal{G}}(n,i)$ for every $d\leq
i\leq d+j$ and with probability at least $1-Cj/d$, $G_{d}\subseteq
G_{d+1}\subseteq\cdots\subseteq G_{d+j}$. Then Theorem 7 follows by taking
$j=\lfloor\epsilon_{n}d\rfloor$.

We prove by induction on $j$. The base case $j=1$ follows directly by Theorem
18 with our choice of $\alpha$. Suppose the statement holds for some $1\leq
j<\lfloor\epsilon_{n}d\rfloor$. Let $\pi_{j}$ be the joint probability
distribution of such a coupling $(G_{d},\ldots,G_{d+j})$. Again, by Theorem
18, there is a coupling $(G_{d+j},G_{d+j+1})$ where
$G_{d+j}\sim{\mathcal{G}}(n,d+j)$, $G_{d+j+1}\sim{\mathcal{G}}(n,d+j+1)$ and
with probability at least $1-C/d$, $G_{d+j}\subseteq G_{d+j+1}$. Let $\pi$
denote the joint probability of this coupling $(G_{d+j},G_{d+j+1})$. We
construct a coupling $(G_{d},\ldots,G_{d+j+1})$ by first sample
$(G_{d},\ldots,G_{d+j})$ according to the distribution $\pi_{j}$, and then
sample $G_{d+j+1}$ according to the conditional probability $\pi(G_{d+j+1}\mid
G_{d+j})$. The resulting coupling satisfies the required marginal distribution
conditions. Moreover, the probability that either
$G_{d}\subseteq\cdots\subseteq G_{d+j}$ fails or $G_{d+j}\subseteq G_{d+j+1}$
fails is at most $Cj/d+C/d=C(j+1)/d$ by the union bound. The assertion follows
by induction.

##  5 Future research

As mentioned in Remark 15, the error $d^{3}/n$ in Corollary 14 can be improved
to $d^{2}/n$ if we apply [4, Theorem 6] and go through more involved
calculations. Another approach is to improve the error in Theorem 13, which is
of independent interest and can lead to improvements of many other existing
results on subgraphs of ${\mathcal{G}}(n,d)$.

We solved Conjecture 5 for a certain range of $d_{1}$, by simultaneously
coupling a sequence of random regular graphs. This is certainly not necessary,
and is the cause of the restrictions on $d_{1}$ in Theorem 6. A more direct
approach would be to prove concentration of the number of $k$-factors in
${\mathcal{G}}(n,d)$. This would significantly relax the restrictions on
$d_{1}$, and itself has independent interest.

## References

  * [1] Béla Bollobás.  A probabilistic proof of an asymptotic formula for the number of labelled regular graphs.  European Journal of Combinatorics, 1(4):311–316, 1980.
  * [2] Andrzej Dudek, Alan Frieze, Andrzej Ruciński, and Matas Šileikis.  Embedding the erdős–rényi hypergraph into the random regular hypergraph and hamiltonicity.  Journal of Combinatorial Theory, Series B, 122:719–740, 2017.
  * [3] Pu Gao.  Kim-vu’s sandwich conjecture is true for all $d=\omega(\log^{7}n)$.  arXiv preprint arXiv:2011.09449, 2020.
  * [4] Pu Gao.  Triangles and subgraph probabilities in random regular graphs.  arXiv preprint arXiv:2012.01492, 2020.
  * [5] Pu Gao, Mikhail Isaev, and Brendan D McKay.  Sandwiching random regular graphs between binomial random graphs.  In Proceedings of the Fourteenth Annual ACM-SIAM Symposium on Discrete Algorithms, pages 690–701. SIAM, 2020.
  * [6] Pu Gao and Yuval Ohapkin.  Subgraph probability of random graphs with specified degrees and applications to chromatic number and connectivity.  arXiv preprint arXiv:2007.02216, 2020.
  * [7] Pu Gao and Cristiane M Sato.  A transition of limiting distributions of large matchings in random graphs.  Journal of Combinatorial Theory, Series B, 116:57–86, 2016.
  * [8] Zhicheng Gao and Nicholas C Wormald.  Distribution of subgraphs of random regular graphs.  Random Structures & Algorithms, 32(1):38–48, 2008.
  * [9] Hans Garmo.  The asymptotic distribution of long cycles in random regular graphs.  Random Structures & Algorithms, 15(1):43–92, 1999.
  * [10] Catherine Greenhill, Mikhail Isaev, and Gary Liang.  Spanning trees in random regular uniform hypergraphs.  arXiv preprint arXiv:2005.07350, 2020.
  * [11] Svante Janson.  The numbers of spanning trees, hamilton cycles and perfect matchings in a random graph.  Combinatorics, Probability and Computing, 3(1):97–126, 1994.
  * [12] Svante Janson.  Orthogonal decompositions and functional limit theorems for random graph statistics, volume 534.  American Mathematical Soc., 1994.
  * [13] Svante Janson.  Random regular graphs: asymptotic distributions and contiguity.  Combinatorics, Probability and Computing, 4(4):369–405, 1995.
  * [14] Jeong Han Kim and Van H Vu.  Sandwiching random graphs: universality between random graph models.  Advances in Mathematics, 188(2):444–469, 2004.
  * [15] Tereza Klimošová, Christian Reiher, Andrzej Ruciński, and Matas Šileikis.  Sandwiching biregular random graphs.  arXiv preprint arXiv:2010.15751, 2020.
  * [16] Brendan D McKay.  Subgraphs of random graphs with specified degrees.  1981\.
  * [17] Brendan D McKay.  Asymptotics for symmetric 0-1 matrices with prescribed row sums.  Ars Combin, 19:15–25, 1985.
  * [18] Brendan D McKay, Nicholas C Wormald, and Beata Wysocka.  Short cycles in random regular graphs.  the electronic journal of combinatorics, pages R66–R66, 2004.
  * [19] Robert W. Robinson and Nicholas C. Wormald.  Almost all cubic graphs are hamiltonian.  Random Structures & Algorithms, 3(2):117–125, 1992.
  * [20] Robert W. Robinson and Nicholas C. Wormald.  Almost all regular graphs are hamiltonian.  Random Structures & Algorithms, 5(2):363–374, 1994.
  * [21] Andrzej Ruciński.  When are small subgraphs of a random graph normally distributed?  Probability Theory and Related Fields, 78(1):1–10, 1988.
  * [22] Nicholas C Wormald.  The asymptotic distribution of short cycles in random regular graphs.  Journal of Combinatorial Theory, Series B, 31(2):168–182, 1981\.
  * [23] Nicholas C Wormald et al.  Models of random regular graphs.  London Mathematical Society Lecture Note Series, pages 239–298, 1999.

## Appendix

  * •

Justify (10)

  * •

Justify (9)

  * •

Justify (17)

We verify

$\sqrt{\frac{n}{e\bar{\delta}(\bar{\alpha}n/2)(n-\bar{\alpha}n)}}\rho_{2}(n,d,\bar{\alpha})=\left(1+\frac{1}{6d^{3}}+O(\xi)\right)\rho_{1}(n,d)^{2}(\bar{\alpha}n/e)^{\bar{\alpha}n/2}.$
(32)

Recall that $\xi=d^{-4}+\frac{d^{3}}{n}+\sqrt{\frac{d}{n}}\log^{6}n$ and
$\bar{\alpha}=1/d$ and

$\rho_{1}(n,d)=\left(\frac{e}{n}\right)^{n/2}\left(\frac{d-1}{d}\right)^{(\frac{d-1}{2})n}d^{\frac{n}{2}}\exp\left(\frac{1}{4}\right),$

and

$\rho_{2}(n,d,\bar{\alpha})=\left(\frac{e}{n}\right)^{(1-\frac{\alpha}{2})n}\left(\frac{d-2+\alpha}{d}\right)^{(\frac{d-2+\alpha}{2})n}d^{\frac{\alpha}{2}n}\left(d-1\right)^{(1-\alpha)n}\exp\left(\phi(d,\bar{\alpha})+O(n^{-1})\right).$

it is easy to check that all exponential terms cancel exactly from both sides
of (32). By Corollary 14 with $\alpha=\bar{\alpha}$ (See Maple calculations
and expansions below),

$\phi(d,\bar{\alpha})=\frac{4d^{2}-10d+5}{4(d-1)^{2}}=1-\frac{1}{2d}-\frac{3}{4d^{2}}-\frac{1}{d^{3}}+O(d^{-4}).$

The polynomially bounded term on the left hand side of (32) is

$\displaystyle\sqrt{\frac{n}{e\bar{\delta}(\bar{\alpha}n/2)(n-\bar{\alpha}n)}}\exp\left(\phi(d,\bar{\alpha})\right)=\sqrt{\frac{d-1}{e(d-2)}}\exp\left(1-\frac{1}{2d}-\frac{3}{4d^{2}}-\frac{1}{d^{3}}+O(d^{-4})\right)$
$\displaystyle\hskip
14.22636pt=\exp(-1/2)\exp\left(\frac{1}{2d}+\frac{3}{4d^{2}}+\frac{7}{6d^{3}}+O(d^{-4})\right)\exp\left(1-\frac{1}{2d}-\frac{3}{4d^{2}}-\frac{1}{d^{3}}+O(d^{-4})\right)$
$\displaystyle\hskip
14.22636pt=\exp\left(\frac{1}{2}+\frac{1}{6d^{3}}+O(\xi)\right).$

See Maple expansion of $\sqrt{(d-1)/(d-2)}$ below where $x=1/d$:

The polynomially bounded term on the right hand side of (32) is

$\left(1+\frac{1}{6d^{3}}+O(\xi)\right)\exp(1/2)=\exp\left(\frac{1}{2}+\frac{1}{6d^{3}}+O(\xi)\right).$