# High precision measurement of the 6He half-life

M. Kanafani  Normandie Université, ENSICAEN, UNICAEN, CNRS/IN2P3, LPC Caen,
14000 Caen, France  X. Fléchard<EMAIL_ADDRESS>Normandie
Université, ENSICAEN, UNICAEN, CNRS/IN2P3, LPC Caen, 14000 Caen, France  O.
Naviliat-Cuncic  Normandie Université, ENSICAEN, UNICAEN, CNRS/IN2P3, LPC
Caen, 14000 Caen, France  National Superconducting Cyclotron Laboratory and
Department of Physics and Astronomy, Michigan State University, East Lansing
48824 MI, USA  G.D. Chung  Normandie Université, ENSICAEN, UNICAEN,
CNRS/IN2P3, LPC Caen, 14000 Caen, France  S. Leblond  Université Paris-Saclay,
CEA, List, Laboratoire National Henri Becquerel (LNE-LNHB), F-91120 Palaiseau,
France  E. Liénard  Normandie Université, ENSICAEN, UNICAEN, CNRS/IN2P3, LPC
Caen, 14000 Caen, France  X. Mougeot  Université Paris-Saclay, CEA, List,
Laboratoire National Henri Becquerel (LNE-LNHB), F-91120 Palaiseau, France  G.
Quéméner  Normandie Université, ENSICAEN, UNICAEN, CNRS/IN2P3, LPC Caen, 14000
Caen, France  A. Simancas Di Filippo  Normandie Université, ENSICAEN, UNICAEN,
CNRS/IN2P3, LPC Caen, 14000 Caen, France  Deutsches Elektronen-Synchrotron
DESY, Notkestraße 85, 22607 Hamburg, Germany  J-C. Thomas  GANIL, CEA/DRF-
CNRS/IN2P3, Bd Henri Becquerel, 14076 Caen, France

###### Abstract

The half-life of 6He has been measured using a low energy radioactive beam
implanted in a YAP scintillator and recording decay events in a 4$\pi$
geometry. Events were time-stamped with a digital data acquisition system
enabling a reliable control of dead-time effects and detector gain variations.
The result, $T_{1/2}=(807.25\pm 0.16_{\rm stat}\pm 0.11_{\rm sys}$) ms,
provides the most precise value obtained so far and is consistent with the
only previous measurement having a precision smaller than 0.1%. This resolves
the longstanding discrepancy previously observed between two sets of
measurements.

††preprint: He6_halflife_pp

##  I Introduction

The study of the 6He beta decay has been instrumental in establishing the V-A
character of the weak interaction Johnson  _et al._ (1963); Vise and Rustad
(1963). The simplicity of the decay has regained considerable attention in the
past decades as a sensitive window to probe the electroweak standard model in
the sector involving the lightest quarks. Experiments using trapped 6He+ ions
Fléchard  _et al._ (2011) or atoms Müller  _et al._ (2022) have been performed
to address the beta-neutrino angular correlation and ancillary measurements
have been carried out to study atomic effects following beta decay Couratin
_et al._ (2012); Hong  _et al._ (2017). New techniques have been proposed to
determine the Fierz interference term from the beta energy spectrum Naviliat-
Cuncic (2016); Huyan  _et al._ (2018); Byron and Helium6 CRES Collaboration
(2020) and selected observables of the transition have been the object of
detailed theoretical studies to assess the sensitivity of future experiments
to new physics Glick-Magid  _et al._ (2021).

One of the properties that makes the 6He decay attractive is the fact that the
transition to the 6Li ground-state is pure Gamow-Teller (GT). As a result, the
correlation terms in the decay rate distribution or the shape of the beta
energy spectrum are, to first order, independent of the dominant nuclear
matrix element, noted $c$ in Ref. Calaprice (1975). However, due to induced
form factors in the hadronic weak currents, this matrix element enters recoil
order corrections in ratios with the weak magnetism form factor Calaprice
(1975). These corrections are important for instance in the extraction of the
Fierz interference term from measurements of the shape of the beta-energy
spectrum.

The value of $c$ extracted from Ref. Knecht  _et al._ (2012) has a relative
uncertainty of $3.6\times 10^{-4}$ so that taking only this value would result
in a negligible impact for the extraction of the Fierz term from a measurement
of the spectrum. However, as pointed out in Ref. Knecht  _et al._ (2012), the
prior internal consistency of the 6He half-lives, having an uncertainty
smaller than 1% is very poor. A fit to the five most precise values prior to
2012 (Fig. 14) gives $T_{1/2}=800.6\pm 0.8$ with $\chi^{2}/\nu=6.95$ with
$\nu=4$. Resolving this discrepancy requires a new measurement with comparable
precision than the one reported in Refs. Knecht  _et al._ (2012) and
preferably with a technique having different sensitivity to the main sources
of systematic effects.

We report here a high precision measurement of the 6He half-life using a
digital data acquisition system which provides a time stamp when recording the
energy of the beta particle. The setup includes an 241Am source which is
constantly monitored in order to control gain variations.

##  II Experimental setup

###  II.1 Beam transport and time sequence

The experiment was performed at the Grand Accélérateur National d’Ions Lourds
(GANIL), Caen. The 6He+ ions were produced by the SPIRAL target-ECR ion source
system and were guided at 25 keV to the low energy beam line LIRAT after mass
separation by a dipole magnet. The beam was chopped with a fast electrostatic
deflector located upstream from the last dipole of the beam line. The
electrostatic deflector was controlled with a Standford Research Systems DG645
pulse generator which also synchronized the motion of the detection system
(see below) and the data acquisition. The vacuum chamber (Fig. 1) is split in
two sections with independent pumping systems. In order to reduce
contamination by 6He atoms, a stainless steel plate with a Ø35 mm hole
separates the two sections. The beam was transported through two Ø6 mm
collimators located at the entrance and exit of the first section of the
chamber [Fig. 1 panel (a)]. During beam tuning, a movable Si detector was
inserted at the entrance of the detection section to measure the incoming beam
intensity.

The detection setup was designed for precision measurements of both, beta
energy spectra and half-lives of specific radioactive species. A $4\pi$
calorimetry coverage is achieved using two detectors. One of the detectors is
fixed and aligned along the beam line whereas the other is mounted on a fast
actuator controlled by the DG645 pulse generator. The 6He+ beam is implanted
on the front face of the fixed detector at a mean depth of about 130 nm as
determined by TRIM, while the moving detector is located out from the beam
[Fig. 1 (a)]. In this configuration, the beam passes across a third Ø4 mm
collimator located a few mm from the fixed detector surface.

The time sequence consists of an implantation interval of typically $T_{\rm
imp}=2.5$ s, followed by a 1 s waiting interval during which the movable
detector is brought in contact with the fixed detector, and finally a
measuring interval of duration $T_{\rm acq}$. The movable detector remains in
the measuring position [Fig. 1 (b)] during data taking before being lifted up
to start over a new cycle.

Figure 1: Sectional view of the vacuum chamber for the implantation (a) and
data taking (b) configurations. The labels on panel (a) are: 1 and 2– the two
Ø6 mm collimators in the first section of the chamber; 3– movable Si detector;
4 and 5– the moving detector and its mechanical guide; 6– the third Ø4 mm
collimator; 7– the fixed detector. The green arrow indicates the 6He+ beam. On
panel (b), label 8 indicates the implantation region.

This closed geometry with two detectors ensures the full collection of beta
particles emitted by the implanted ions and prevents any partial energy loss
due to backscattering. The duration of the data taking sequence was adjusted
to be long enough such as to obtain a precise measurement of the ambient
background. Here after, the fixed and movable detectors will be referred to as
“Det1” and “Det2” respectively. The transit times of Det2 were measured to be
0.8 s on the way down and 1.7 s on the way up. The loss of activity between
the end of the implantation interval and the beginning of the measurement was
thus limited to 50%. The total transit time of 2.5 s for a cycle duration of
about 20 s does not lead to a significant loss of statistics. The motion of
the detector was intensively tested for several days with 20 seconds duration
cycles prior the experiment. These tests have shown no degradation of the
detector signals. The fair contact between the two detector faces while in
data taking position was also checked before and after the experiment by
pinching 20 $\mu$m thickness sheets at different places of the contact
surface.

###  II.2 Beta particle detectors

The two detectors, Det1 and Det2, are identical. Each of them is composed of a
cylindrical Ø$30\times 30$ mm2 YAlO3 Ce doped inorganic scintillator (YAP)
surrounded by a EJ-204 plastic scintillator (PVT) with external diameter of 40
mm (Fig. 2). The scintillators are mounted in a phoswich configuration in
which the YAP and the PVT are readout by a single R7723 photo-multiplier tube
(PMT) from Hamamatsu. The aluminum detector housing serves in part as a shield
and ensures also some pressure to produce a fair optical contact between the
PMT, the optical coupling grease and the scintillators. Two circular hollows
at the front of the housing allow the insertion of calibration sources.

Figure 2: Left panel: Close sectional view of a detector with the two
scintillators coupled to a single PMT in a phoswich configuration and the
location of the 241Am source. Right panel: Typical signals from the PVT and
YAP along with the three charge integration intervals.

A 5 kBq 241Am source was permanently mounted on each detector to provide a
constant calibration reference using the 59.54 keV photons. These interact
mostly in the YAP volume whereas the $\sim$5 MeV alpha particles are stopped
in the PVT which also serves as a veto to reject background events. The
detection of the alpha particles by the PVT is used to monitor PMT gain
variations. A $225~{}\mu$m layer of Tyvek with a Ø5 mm hole in front of the
241Am source was wrapped around the PVT to improve the light collection. The
faces of both detectors which get into contact are free of any reflector or
dead layer. As a result, when scintillation light is generated in one of the
detectors, about 80% of the light is collected by its own PMT while the other
20% is collected by the PMT of the other detector.

The choice of YAP crystals resulted from a compromise between a fast response,
an acceptable energy resolution as well as differential linearity over the
largest possible energy range Moszyński  _et al._ (1998); Mengesha  _et al._
(1998). The YAP time decay constant of 27 ns is short enough to limit pile-up
contributions but significantly longer than the PVT decay time of 2 ns. This
leads to a clear discrimination between the PVT and YAP signals by pulse-shape
analysis. Typical signals from the PVT and from the YAP scintillators are
shown on the right panel of Fig. 2. Moreover, the high photon detection
efficiency of YAP allows the use of the photopeaks from gamma sources for off-
line energy calibrations. The nominal PMT polarization voltages used for Det1
and Det2 were $-1620$ V and $-1420$ V respectively. These values were chosen
in order to ensure a linear response of the detectors over the full beta
energy range (up to 3.5 MeV) and also to match their relative gain.

The response functions of the detectors were studied in detail prior to
installing the detectors inside the vacuum chamber, using 241Am, 22Na, 137Cs
and 60Co sources, as well as $\gamma$ rays from 208Tl. Small PMT gain
variations between the measurements with different sources were corrected
using the reference provided by the alpha particles interacting in the PVT.
Deviations from linearity were found to be smaller than 4 keV over the 60
keV–2.6 MeV range covered by the calibration sources. However, these
deviations are most likely due to the fit of the photopeak which neglects the
presence of the overlapping Compton or multi-Compton distributions. The energy
resolution was 6.7% at 662 keV and 5.5% at 1332 keV.

Before and after the experiment, data were recorded again with the detectors
under the experimental conditions and using 241Am, 22Na, 137Cs and 60Co
sources to check the response function of the detectors.

###  II.3 Data acquisition system

The analog signals from the two PMTs and a logic pulse from the DG645
generator were directly sent to three input channels of the digital data
acquisition system FASTER FAS (p3fr). Each channel digitizes the signal at a
rate of 500 MS/s with each sample converted over 12 bits. The triggers of the
channels are independent based on individual thresholds. The digitized samples
are processed in real time by Field-Programmable Gate Arrays using predefined
algorithms adapted to the measurements to be performed FAS (p3fr). All samples
are time-stamped with a 2 ns resolution allowing on-line and off-line
correlations within user-defined time windows. For the two PMTs, the selected
algorithm provided the charge integration of the signals within four different
time windows relative to the trigger time. The first window was set between
$-30$ and $-6$ ns to obtain the charge $Q_{\rm baseline}$ which provides an
estimate of the baseline level (Fig. 2 right panel). The second and third
windows were set between $-6$ and 30 ns and between $-6$ and 300 ns
respectively. They enable to cover the PVT fast signals (charge $Q_{\rm
fast}$) and the slower signals from the YAP (charge $Q_{\rm tot}$). The ratio
$Q_{\rm fast}/Q_{\rm tot}$ serves for pulse-shape analysis and signal
discrimination. A fourth window between $-30$ ns and 1 $\mu$s was also set to
study after-pulses and pile-up effects. With these settings, the intrinsic
dead-time of the acquisition system was 1030 ns. The logic pulse from the
DG645 generator was used to provide a time reference for each measurement
cycle. The baseline of the three signals was continuously monitored for each
channel and corrected for low-frequency variations (below 160 kHz) by the
FASTER baseline restoration algorithm. The trigger threshold of the two PMT
channels was set to 5 mV, which corresponds to a minimum deposited energy of
about 4 keV. The stability and accuracy of the internal clock of the FASTER
system is smaller than 1 ppm and the associated systematic uncertainty on the
half-life measurement is thus negligible relative to other corrections.

##  III DATA ANALYSIS

In order to study systematic and background effects, five different
experimental conditions were adopted resulting in five statistically
independent data sets. The three parameters which were changed among the sets
were the cycle duration, $T_{\rm acq}$, the initial activity and the bias
voltages of the PMTs. A measurement of the ambient background with beam ON but
without any implantation was also performed using a 0.4 mm thick aluminum disk
fully covering the third collimator. For these runs, the incident beam
intensity was adjusted so as to match the conditions of the short and long
cycles runs at nominal PMT voltages, defining the conditions of sets (4) and
(5) respectively. The different sets and the corresponding number of cycles
are summarized in Table 1.

Table 1:  The five experimental conditions of the recorded data. $T_{\rm acq}$ is the duration of the measuring window; $\langle N_{0}\rangle$ is the average number of detected 6He decays per cycle with an energy threshold of 300 keV; $\Delta V_{\rm PMT}$ indicates the difference in the PMT bias relative to the nominal voltage. The last column lists the number of recorded cycles.  Set | Set name |  $T_{\rm acq}$ (s) | $\langle N_{0}\rangle$ |  $\Delta V_{\rm PMT}$ (V) | Cycles
---|---|---|---|---|---
(1) | Short cycles | 11 | 6130 | 0 | 3050
(2) | Long cycles | 26 | 20900 | 0 | 620
(3) | Lower bias/Short | 11 | 7140 | $-50$ | 1930
(4) | Shutter/Short | 11 | $-$ | 0 | 200
(5) | Shutter/Long | 26 | $-$ | 0 | 115

For sets (1), (2) and (3), runs were devoted at the beginning of each
measurement to record signals waveforms which provided evidence of the
presence of after-pulses in the range between 400 ns and 1 $\mu$s after the
trigger. This prevented the use of the fourth charge integration window for
the extraction of the 6He half-life.

###  III.1 Baseline and gain corrections

A precise absolute energy calibration is not required for a precise half-life
measurement. However, it is crucial to control baseline and gain variations
within the decay cycle which could result in a change of the effective energy
threshold imposed to the data.

The method described below is based on monitoring the time dependence of the
59.54 keV line from the 241Am source. The calibrations of Det1 and Det2 were
performed independently for each set of measurements. The method is
illustrated here with data from set (1) and for Det1. First, a 2D-plot of the
ratio $Q_{\rm fast}/Q_{\rm tot}$ versus $Q_{\rm tot}$ is built which enables
particle identification through pulse-shape discrimination. Events associated
with scintillation light generated by the YAP only are then selected by
setting a window over $Q_{\rm fast}/Q_{\rm tot}$ between 0.5 and 0.62 (Fig.
3).

Figure 3: 2D-plot of the the ratio $Q_{\rm fast}/Q_{\rm tot}$ versus the total
charge for signals from Det1. The contributions of the different particles are
indicated. The horizontal red lines define the window of events due to an
interaction in the YAP volume only. The plot was obtained for a one hour run
under conditions of set (1).

The cycles were then sorted as a function of the signal to background ratio
(SBR) defined as the number of 6He decay events divided by the average
constant background measured within a single cycle. The SBR was obtained for
each cycle using a fit of the time distribution with an exponential decay plus
a constant. The SBR distribution for set (1) is shown in Fig. 4. Five groups
of SBR values containing comparable number of cycles were then defined. Cycles
with SBR smaller than 5.0 were definitely discarded from the analysis as their
very low statistics can lead to fitting convergence problems in the half-life
extraction. The small cluster at SBR $\sim 25$ in Fig. 4 was also excluded
from the baseline and gain correction model because of its low statistics, but
the corresponding cycles were included in the data used for the final half-
life estimate.

Figure 4: Distribution of SBR for set (1). The groups for different values of
SBR, labelled from 1 to 5, are indicated by the vertical red lines. Figure 5:
Events collected during a one hour long run under conditions of set (1) as a
function of the time within a cycle. The vertical red lines indicate the
separation between the time windows selected for the gain variation
correction. The beginning of the beam-on period is at $t=0$ and the
measurement period starts at $t=3.5$ s. The blue histogram corresponds to the
total number of events, including those from the 241Am source. The green
histogram are those events remaining after the $Q_{\rm fast}/Q_{\rm tot}$
charge selection and an energy threshold of 100 keV. One bin corresponds to 50
ms.

For each group of cycles, events were then sorted according to their detection
time within the cycle, using successive windows of 1 s duration (Fig. 5). For
each set of events corresponding to a given time window and a given SBR
window, the charge distribution $Q_{\rm tot}$ of the 59.54 keV photopeak and
the $Q_{\rm baseline}$ distributions were fitted using Gaussian functions
[Fig. 6 (a) and (b)]. The mean values ($M_{\rm tot}$ and $M_{\rm baseline}$)
obtained from these fits are plotted in Fig. 6 (c) and (d) as a function of
time for the five different windows of SBR. A significant gain variation,
close to 2%, is observed in Fig. 6 (c) between the beginning and the end of
the measurement period. A baseline variation is also visible but has a smaller
impact. For each SBR window, the data of Fig. 6 (c) and (d) were then fitted
with functions of the form $P_{0}+P_{1}e^{-t/\tau_{G}}$ where $\tau_{G}$ is a
common parameter for the five different SBR windows.

Figure 6: (a) and (b): Close-up to the $Q_{\rm tot}$ and $Q_{\rm baseline}$
distributions associated to a SBR window and a time window. The red lines show
the fits of the 59.54 keV photopeak (a) and of the baseline (b). (c) and (d):
The data points show the mean values resulting from the fits as a function of
time for the 59.54 keV photopeak (c) and the baseline (d). The black, red,
green, blue and yellow colors correspond respectively to the SBR windows 1 to
5. The lines are the exponential fit functions obtained for each SBR set.

The gain and baseline correction models were completed by extracting
parameters $P_{0}$ and $P_{1}$ from fits of the mean values of the charge
distributions as a function of time for each SBR window. For the gain
correction, the parameters are illustrated in Fig. 7 (a) and (b). The
parameter $P_{1}$ is here the most relevant one since it causes a time
dependence of the gain and of the baseline. It was found to depend linearly on
the SBR value for both the gain and baseline distributions. For the gain
correction model, $P_{0}$ was extracted using an exponential function whereas
it was taken as constant for the baseline correction model.

Figure 7: (a) and (b): Parameters $P_{1}$ and $P_{0}$ for the gain correction
model for the five values of SBR; (c) Results from fits of the 59.54 keV
photopeak as a function of time without (magenta) and with (blue) the gain and
baseline corrections.

Using these models, Det1 and Det2 were independently calibrated for each time
window within a cycle and using the parameters associated with the SBR value
of the cycle. Calibrated data from both detectors were then summed up for each
event and scaled by a final calibration coefficient of 0.82, evaluated run by
run, which accounts for the light cross-talk between the detectors for photons
of 60 keV. Fig. 7 (c) shows an example of the values from the fits of the
59.54 keV photopeak, as a function of time, with (blue) and without (magenta)
applying the gain and baseline corrections.

###  III.2 Background sources

The background rate was assumed to be constant during a cycle. The measuring
cycle was divided in two halves and the energy spectrum corresponding to decay
events was then obtained by subtracting the spectrum recorded during the
second half from the spectrum recorded during the first half. The resulting
energy spectra obtained for sets (1) and (4) are shown in Fig. 8. The black
histogram corresponds to set (1) and shows two contributions. The first is the
beta energy spectrum extending up to 3.5 MeV and the second is a low energy
distribution with a peak at about 0.1 MeV. The decay time of this distribution
is consistent with the 6He half-life and thus the contribution was attributed
to Bremsstrahlung radiation of electrons from 6He+ decay, for ions implanted
in the third collimator. This was confirmed by the analysis of data from set
(4), shown by the green histogram in Fig. 8. We recall that the data for this
set was obtained with a shutter located inside the third collimator in order
to prevent any implantation on the detector. The theoretical beta energy
spectrum convoluted with the detector response function, is shown by the blue
histogram. The sum of the blue and green distributions is displayed by the red
histogram which fairly reproduces the data from set (1) (insert in Fig. 8).
This Bremsstrahlung background does not affect the measurement since it decays
with the same half-life. Furthermore, by setting an energy threshold above 0.6
MeV, its contribution represents 0.5% of the number of events associated with
the beta particle detection.

Figure 8: Energy spectra for set (1) (black histogram) and (4) (green
histogram) after background subtraction. The blue histogram is the theoretical
beta energy spectrum including the convolution with the detector response
function. The red histogram is the fit to the black spectrum with a
combination from the blue and green spectra. The green and blue distributions
are normalized using the parameters from the fit.

Conversely, the selection of data in the second half of the measurement cycle
suppresses the contribution due to 6He decay. In order to fully subtract this
contribution, a proper normalization factor was used on the spectrum measured
during the first half of the cycle. This normalization factor was determined
assuming that events with an energy deposit larger than 2 MeV were solely due
to 6He decay. The resulting energy spectrum is shown in Fig. 9 (red histogram)
and corresponds to a constant or long half-life background compared to 6He.
The peak at 59.54 keV is from the monitoring 241Am source and the photo-peaks
at 662 keV and 1461 keV are from a 137Cs source of a dose rate meter located
in the experimental area and from ambient 40K respectively. Compared to the
ambient background spectrum recorded before and after the on-line experiment
(black histogram), an excess of events is observed below 0.6 MeV, with a peak
due to 511 keV photons. This contaminant was not clearly identified and its
half-life was found too long to be measured with the cycle length chosen for
the experiment. A lower limit of 330 s (95% CL) could however be extracted
from this background analysis. The average detection rate of this contaminant
contribution was found to be $\sim$14 cps and it can be fully suppressed by
applying an energy threshold of 0.6 MeV. The ambient background contribution
was constant during the experiment.

Figure 9: Energy spectra for ambient background (black histogram) and for
constant or long half-life background recorded on-line for one run of data set
(1) (red histogram). The black histogram was normalized to the red one by
using the 59.54 keV peak from 241Am.

###  III.3 Dead-time correction

To study systematic effects due to dead time, a common recipe Knecht  _et al._
(2012) consists in using analog parallel channels with different non-
extensible dead times in the data acquisition. With precisely time-stamped
events recorded by the digital system, the data were first filtered off-line
by imposing several dead times $\tau_{DT}$ ranging from 1.03 to 7 $\mu$s. In
order to study further possible bias due to background contributions or
imperfections of the gain and baseline corrections, the data were additionally
filtered using energy thresholds ranging from 0.1 to 1.2 MeV.

The probability to miss a decay event at time $t$ within a cycle because of
dead time can be expressed as $\tau_{DT}\cdot r_{T}(t)_{i}$, where $\tau_{DT}$
is the non-extensible dead time and $r_{T}(t)_{i}$ is the total instantaneous
rate of events that could potentially trigger the data acquisition system
during cycle $i$. This rate is given by

$r_{T}(t)_{i}=r_{0i}e^{-t/\tau_{i}}+r_{bi}~{},$ (1)

with $r_{0i}$ the total initial decay rate and $r_{bi}$ the total rate of
background events. First, the parameters $r_{0i}$, $\tau_{i}$ and $r_{bi}$
which enter Eq.(1) were obtained by fitting the rate of detected events for
each cycle, $r_{D}(t)_{i}$, without threshold condition, and using the
correction accounting for losses due to dead time

$r_{D}(t)_{i}=\frac{r_{T}(t)_{i}}{1+\tau_{DT}\cdot r_{T}(t)_{i}}~{}.$ (2)

If only the dead-time effect was present, for a given dead time $\tau_{DT}$,
the occurrence of an event at time $t$ should be weighted by a coefficient

$w(t)_{i}=1+\tau_{DT}\cdot r_{T}(t)_{i}~{}.$ (3)

However, pile-up effects have also to be considered and can be included in the
weighting coefficient.

###  III.4 Pile-up effects

An off-line energy threshold condition, labelled $j$, was also applied to the
data. With this condition, the rate of selected events is then given by

$r_{D}(t)_{ij}=r_{D}(t)_{i}\cdot P_{E}(t)_{j}~{},$ (4)

where $P_{E}(t)_{j}$ is the probability for the measured deposited energy to
be above the energy threshold $E_{j}$. To first order, $P_{E}(t)_{j}$ should
be independent of time but the occurrence of pile-up events during the charge
integration window can lead to a small time dependent correction. When a pile-
up occurs within the window, the recorded charge and its associated energy
will naturally be larger than for each of the events producing the pile-up.
The recorded event has then a larger probability to result in an energy above
the threshold. The pile-up probability is proportional to the instantaneous
rate in the decay cycle. This leads to an excess of events with an energy
above threshold for high rates as compared to low rates and would result in an
underestimation of the half-life if not properly corrected.

The instantaneous energy-independent probability for a decay event to pile-up
with another decay event within a cycle $i$ is given by

$P_{\rm pu}(t)_{i}=r_{0i}e^{-t/\tau_{i}}\Delta T~{},$ (5)

where $r_{0i}$ and $\tau_{i}$ are the parameters determined from the fit of
Eq.(1) and $\Delta T=300$ ns is the duration of the charge integration window.

The excess of detected events due to pile-up can then be expressed as

$r_{E}(t)_{ij}\approx P_{\rm
pu}(t)_{i}\cdot(d_{j}-s_{j})r_{0i}e^{-t/\tau_{i}}~{},$ (6)

where $s_{j}$ and $d_{j}$ are respectively the time independent probabilities
for an energy conversion to be above threshold for single events and for pile-
up events involving two signals. The contribution to the half-life of pile-up
events involving three or more signals was estimated to be smaller than
$10^{-7}$ s and hence negligible. On the other hand, the expected rate without
pile-up contribution would be

$r_{R}(t)_{ij}=s_{j}\cdot r_{0i}e^{-t/\tau_{i}}+b_{j}\cdot r_{bi}~{},$ (7)

where $b_{j}$ is the time independent probability to be above threshold for
background events.

The probabilities $s_{j}$ and $b_{j}$ to get a signal above threshold $E_{j}$
for single events were calculated from the deposited energy distribution of
6He decay events and background events, respectively (Figs. 8 and 9). The
expected energy distribution for 6He pile-up events was deduced using the
auto-convolution of the energy distribution obtained for 6He single events.
The mean probability, $d_{j}$, for a pile-up event within 300 ns to be above
threshold was then determined. The fraction of the charge of the second signal
that is lost when overlapping with the end of the integration window was
accounted for in this calculation. Note that the effect of pile-up events
involving the constant ambient background and the contaminant contribution
discussed in section III.B were neglected in Eq.(6) since the associated
correction was estimated to be below 10-6 s.

Finally, the relative rate excess of detected events as a function of time can
then be approximated by

$\alpha_{E}(t)_{ij}=\frac{r_{E}(t)_{ij}}{r_{R}(t)_{ij}}~{},$ (8)

and the weighting coefficient which accounts for both, dead-time and pile-up
effects becomes

$w(t)_{ij}=\frac{1+\tau_{DT}\cdot r_{T}(t)_{i}}{1+\alpha_{E}(t)_{ij}}~{}.$ (9)

For simplicity in the discussion above, a given value of the dead time,
$\tau_{DT}$, was taken in all expressions. Changing this value corresponds to
an additional condition, labelled $k$ here below, so that the weights in
Eq.(9) become $w(t)_{ijk}$.

###  III.5 Fit procedure

After applying a given dead time and energy threshold, the weighted events
from all cycles within a measurement set were summed-up and binned. The
weighted number of counts in a bin is then $n(t)_{jk}=\sum_{i}{w(t)_{ijk}}$
and the variance is $\sigma^{2}_{jk}=\sum_{i}{w^{2}(t)_{ijk}}$, where the sums
run over the cycles and also over all events within each cycle. The corrected
data were finally fitted assuming a constant background. For a given dead time
and energy threshold, the fit function can be expressed as

$f(t)_{jk}=A_{jk}\cdot e^{-t/\tau_{jk}}+B_{jk}~{},$ (10)

where $A_{jk}$ is the initial number of decay counts, $\tau_{jk}$ is the
estimate of the decay lifetime and $B_{jk}$ is the constant background level.

To summarize, the fit procedure for each set involves three steps. First, the
rate of detected events without energy condition was fitted for each cycle
with Eq.(2) with $r_{0i}$, $\tau_{i}$ and $r_{bi}$ as free parameters, to be
used in the dead-time and pile-up corrections. Next, the dead time and energy
threshold corrections were estimated and the occurrence of each event was
weighted using Eq.(9). Finally, weighted events from a measurement set were
summed-up and the resulting binned distribution was fitted with Eq.(10), with
$A_{jk}$, $\tau_{jk}$ and $B_{jk}$ as free parameters. Following this
procedure, all the data from a measurement set results in a single histogram.
For all fits, the parameter estimates were obtained using the log-likelihood
function. The nominal bin size adopted in the histograms was 50 ms.

Figure 10 shows the fit of data set (2) with a dead time of 3 $\mu$s and an
energy threshold of 600 keV. The fit of the standard residuals distribution is
consistent with a normal distribution with $\mu=-0.027(37)$ and
$\sigma=1.01(4)$.

Figure 10: Upper panel: experimental decay spectrum from set (2) for a non-
extensible dead time of 3 $\mu$s and an energy threshold of 600 keV (black)
along with its fit function (red). Lower panel: standard residuals. The
standard residuals distribution and its fit by a Gaussian are shown in the
insert of the upper panel.

In Table 1, the sets from which the half-life can be extracted are (1), (2)
and (3). The resulting half-life values are shown in Fig. 11 as a function of
the imposed dead time and for three energy thresholds. The values show no
significant effect due to dead time except for set (2) and
$\tau_{DT}=1~{}\mu$s, which is slightly lower than the others. This
discrepancy may be due to the contribution of after-pulses in the $1~{}\mu$s -
$2~{}\mu$s range. In the following, a conservative value $\tau_{DT}=3~{}\mu$s
was adopted. The systematic uncertainty on the dead-time correction was
obtained from the 2 ns accuracy on the imposed dead time and is smaller than
10-5 s for all sets. The results from the fits of the histograms for the three
sets are summarized in Table 2. The central values include the systematic
corrections due to dead time and pile-up but the uncertainties are only
statistical.

Figure 11: Mean values of the half-life estimates obtained for set (1) (upper panel), (2) (middle panel) and (3) (lower panel) as a function of the non-extensible dead time. The black squares, red dots and blue triangles correspond respectively to off-line energy thresholds of 300, 600 and 900 keV. The statistical error bars are only shown for the threshold of 600 keV. Table 2: Values of the 6He half-life obtained from the fits of histograms for the three data sets along with their associated $p$-value. | Set (1) | Set (2) | Set (3)
---|---|---|---
$T_{1/2}$ [ms] | 807.42(25) | 807.16(26) | 807.10(35)
$p$-value | 0.70 | 0.83 | 0.25

###  III.6 Effect of baseline and gain corrections

For illustration, Fig. 12 shows, for the three sets, the evolution of the mean
values of the half-life as a function of the energy threshold, before (black
circles) and after (red squares) applying the gain correction described in
Sec. III.1. Whereas the non-corrected data display a trend which increases
with the threshold between 300 and 1200 keV, the corrected data lead to mean
values which are independent of the threshold. The results obtained without
gain correction for an energy threshold of 200 keV are out of the vertical
range of the figure, at about 20 ms higher than with gain correction. This is
due to the strong negative slope of the energy spectrum at 200 keV, right from
the Bremsstrahlung peak shown in Fig. 8. The systematic uncertainty associated
to the gain correction procedure is shown by the gray areas. It was estimated
by accounting for the uncertainties on the parameters $P_{0}$ and $P_{1}$ of
the model and taking values which maximize or minimize the gain correction
amplitude. As expected, this uncertainty also increases with the energy
threshold between 300 and 1200 keV and is very large for a threshold of 200
keV. The distribution of the corrected values as a function of the threshold
is consistent with a constant when accounting for both statistical and
systematic uncertainties. In order to minimize the uncertainty on the gain
correction, one should normally favor a threshold in the range 300-400 keV.
However, because of the presence of background contributing up to 550 keV
whose half-life is unknown (Fig. 9) a threshold of 600 keV was finally chosen.
The effect of the baseline correction was studied in a similar way and was
found to be one order of magnitude smaller than for the gain correction. The
$\chi^{2}$ and residuals distributions obtained for all sets where also
studied with and without gain corrections. In both cases, they showed no
deviation from the expected statistical fluctuations, which indicates that
analyses of $\chi^{2}$ or $p$-values do not provide a proper diagnostic to
detect such systematic effects.

Figure 12: Mean values of the half-life obtained for set (1) (upper panel),
(2) (middle panel) and (3) (lower panel) as a function of the energy
threshold. The black circles (resp. red squares) are the values obtained
without (resp. with) the gain correction. The error bars are statistical. The
horizontal blue line shows the value with correction for a threshold of 600
keV. The gray area indicates the systematic uncertainty associated to the gain
correction and the light gray one the combined systematic and statistic
uncertainty (see text for details).

###  III.7 Diffusion of 6He

The possible rapid diffusion of 6He atoms out of the detector bulk was
considered. Diffusion coefficients specific to helium implanted in a YAP
crystal could not be found in the literature. However, experimental data are
available for a number of mineral compounds Trocellier  _et al._ (2014).
Diffusion coefficients at room temperature are all smaller than 10-26 m2s-1.
Using Fick’s second law, an implantation depth of 100 nm leads to effusion
time constants larger than 1011s. Such time constants are far too large to
have any significant effect on the 6He half-life measurement.

Figure 13: Half-life estimate obtained for the measurements sets (1), (2) and
(3) (in gray, blue and red, respectively) as a function of the lower bound of
the fitting range, $t_{min}$. The colored area indicate the range covered by
the statistical uncertainties given at one standard deviation. The horizontal
black line is the central value of the final result.

###  III.8 Additional cross-checks

The bin width of the histograms was changed from the nominal 50 ms to 20 and
100 ms to check the sensitivity to the binning. The resulting values of the
half-life were found to be the same within 1$\times 10^{-5}$ s.

The dependence of the half-life on the lower bound of the fitting range,
$t_{min}$ was also investigated (Fig.13). The analysis of $p$-values spanning
more than 6 half-lives indicates that the observed variations are fully
consistent with statistical fluctuations within the three sets.

In the analysis procedure presented above, the data was corrected and then
fitted with Eq.(10). An additional analysis was also carried out without
applying any dead time and pile-up corrections to the data. Instead, for a
given condition in dead time and energy threshold, the data in each individual
cycle $i$ was fitted with a function accounting for both dead time and pile-up
effects

$f_{ij}(t)=\frac{(A_{ij}\cdot
e^{-t/\tau_{ij}}+B_{ij})[1+\alpha_{E}(t)_{ij}]}{1+\tau_{DT}\cdot
r_{T}(t)_{i}}~{},$ (11)

where $r_{T}(t)_{i}$ and $\alpha_{E}(t)_{ij}$ are the functions in Eqs.(1) and
(8) determined for each cycle and selection and where $A_{ij}$, $\tau_{ij}$
and $B_{ij}$ are the free parameters of the fit. For each measurement set, the
lifetime estimates, $\tau_{ij}$, obtained from the fits were averaged to
obtain the half-life $\tau_{j}$ for that condition. This averaging required a
closer look because the low statistics present in many cycles was found to
induce a correlation between the half-life and the absolute statistical
uncertainty.

In order to test the averaging procedures, Monte Carlo simulations were
performed using typical experimental decay and background rates and a higher
number of decay cycles. The simulations showed that when using the relative
statistical uncertainty instead of the absolute uncertainty as the weighting
factor, the bias due to low statistics in the fitted data was strongly
reduced. Nevertheless, the weighted average was found to be overestimated by
9$\times 10^{-5}$ s for the conditions of sets (1) and (3) and by 1.5$\times
10^{-5}$ s for the condition of set (2). After accounting for the bias due to
fits with low statistics, the final experimental values obtained when fitting
independently each cycle were found to be the same within 1$\times 10^{-5}$ s
to those obtained using the method described in Sec. III.5.

The complete analysis was performed again while using the result for the half-
life obtained in Eq.(12) to fix the value of $\tau_{i}$ in Eqs.(1) and (2).
The results were again found to be the same within 1$\times 10^{-5}$ s.

##  IV RESULT

Table 2 gives the summary of the values obtained from the fits which include
the main systematic corrections, whereas Table 3 gives the size of the
systematic corrections to those values along with their associated
uncertainties. The values from the three sets are statistically consistent,
with differences in the central values smaller than one standard deviation.
The systematic corrections add up to about 1 ms. The largest shift is due to
the gain variation correction and indicates how crucial it is to control this
effect when aiming at a relative precision smaller than $10^{-3}$. The
combined results of the three measurement sets yield the value

$T_{1/2}=(807.25\pm 0.16_{\rm stat}\pm 0.11_{\rm syst})~{}~{}{\rm ms}~{},$
(12)

where the largest systematic uncertainty has been adopted. Figure 14 shows a
comparison between the present result (horizontal lines) and the six
previously measured values having a relative precision smaller than 1%. The
present result is consistent with three of the previous values Wilkinson and
Alburger (1974); Alburger (1982); Knecht  _et al._ (2012) and is at variance
with three others Kline and Zaffarano (1954); Bienlein and Pleasonton (1962);
Barker  _et al._ (1981). Together with the result from Ref. Knecht  _et al._
(2012), which has a similar precision but used a different technique, the
present result confirms the 6He half-life close to 807 ms and strongly
disfavors previous results yielding values below 800 ms.

Table 3: Corrections to the half-life (in ms) associated with the main sources of systematic effects for the three sets of data. The values of the uncertainties are rounded at 0.01 ms. The combined uncertainties of the total correction were obtained by summation in quadrature. Source | Set (1) | Set (2) | Set (3)
---|---|---|---
Gain | 0.75(7) | 0.77(10) | 0.78(6)
Baseline | 0.09(3) | 0.04(2) | 0.05(9)
Pile-up | 0.10(1) | 0.25(1) | 0.11(1)
Binning | $<0.01$ | $<0.01$ | $<0.01$
Total correction | 0.94(7) | 1.06(11) | 0.94(11)
Figure 14: Comparison between the 6He half-life value obtained from the
present work (horizontal lines) and previous measurements having a relative
precision smaller than 1%. The plotted values are from Refs.Kline and
Zaffarano (1954); Bienlein and Pleasonton (1962); Wilkinson and Alburger
(1974); Barker  _et al._ (1981); Alburger (1982); Knecht  _et al._ (2012).

##  V CONCLUSION

This work reported the most precise value of the 6He half-life obtained so
far. The result is consistent with the most recent measurement which supported
two previous values around 807 ms. By recording both the deposited energy and
the time of each event with a digital data acquisition system combined with
the use of a monitoring 241Am source, detector gain variations and dead-time
effects were precisely measured and corrected for. These effects were found to
contribute at a relative level of 10-3 with a resulting relative systematic
uncertainty at the level of 10-4.

##  VI ACKNOWLEDGMENTS

The authors thank the LPC Caen and GANIL staffs for the technical support and
are grateful to A. Singh for her assistance during the running of the
experiment. This project was supported in part by the French Agence Nationale
de la Recherche under grant ANR-20-CE31-0007-01 (bSTILED).

## References

  * Johnson  _et al._ (1963) C. H. Johnson, F. Pleasonton, and T. A. Carlson, Phys. Rev. 132, 1149 (1963).
  * Vise and Rustad (1963) J. B. Vise and B. M. Rustad, Phys. Rev. 132, 2573 (1963).
  * Fléchard  _et al._ (2011) X. Fléchard, P. Velten, E. Liénard, A. Méry, D. Rodríguez, G. Ban, D. Durand, F. Mauger, O. Naviliat-Cuncic, and J. C. Thomas, Journal of Physics G: Nuclear and Particle Physics 38, 055101 (2011).
  * Müller  _et al._ (2022) P. Müller, Y. Bagdasarova, R. Hong, A. Leredde, K. G. Bailey, X. Fléchard, A. García, B. Graner, A. Knecht, O. Naviliat-Cuncic, T. P. O’Connor, M. G. Sternberg, D. W. Storm, H. E. Swanson, F. Wauters, and D. W. Zumwalt, arXiv (2022), 10.48550/ARXIV.2206.00742.
  * Couratin  _et al._ (2012) C. Couratin, P. Velten, X. Fléchard, E. Liénard, G. Ban, A. Cassimi, P. Delahaye, D. Durand, D. Hennecart, F. Mauger, A. Méry, O. Naviliat-Cuncic, Z. Patyk, D. Rodríguez, K. Siegień-Iwaniuk, and J.-C. Thomas, Phys. Rev. Lett. 108, 243201 (2012).
  * Hong  _et al._ (2017) R. Hong, A. Leredde, Y. Bagdasarova, X. Fléchard, A. García, A. Knecht, P. Müller, O. Naviliat-Cuncic, J. Pedersen, E. Smith, M. Sternberg, D. W. Storm, H. E. Swanson, F. Wauters, and D. Zumwalt, Phys. Rev. A 96, 053411 (2017).
  * Naviliat-Cuncic (2016) O. Naviliat-Cuncic, AIP Conf. Proc. 1753, 060001 (2016).
  * Huyan  _et al._ (2018) X. Huyan, O. Naviliat-Cuncic, P. Voytas, S. Chandavar, M. Hughes, K. Minamisono, and S. Paulauskas, Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment 879, 134 (2018).
  * Byron and Helium6 CRES Collaboration (2020) W. Byron and Helium6 CRES Collaboration, in  _APS Division of Nuclear Physics Meeting Abstracts_ , APS Meeting Abstracts, Vol. 2020 (2020) p. MF.002.
  * Glick-Magid  _et al._ (2021) A. Glick-Magid, C. Forssén, D. Gazda, D. Gazit, P. Gysbers, and P. Navrátil, arXiv (2021), 10.48550/ARXIV.2107.10212.
  * Calaprice (1975) F. P. Calaprice, Phys. Rev. C 12, 2016 (1975).
  * Knecht  _et al._ (2012) A. Knecht, R. Hong, D. W. Zumwalt, B. G. Delbridge, A. García, P. Müller, H. E. Swanson, I. S. Towner, S. Utsuno, W. Williams, and C. Wrede, Phys. Rev. C 86, 035506 (2012).
  * Moszyński  _et al._ (1998) M. Moszyński, M. Kapusta, D. Wolski, W. Klamra, and B. Cederwall, Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment 404, 157 (1998).
  * Mengesha  _et al._ (1998) W. Mengesha, T. Taulbee, B. Rooney, and J. Valentine, IEEE Transactions on Nuclear Science 45, 456 (1998).
  * FAS (p3fr) (https://faster.in2p3.fr/).
  * Trocellier  _et al._ (2014) P. Trocellier, S. Agarwal, and S. Miro, Journal of Nuclear Materials 445, 128 (2014).
  * Wilkinson and Alburger (1974) D. H. Wilkinson and D. E. Alburger, Phys. Rev. C 10, 1993 (1974).
  * Alburger (1982) D. E. Alburger, Phys. Rev. C 26, 252 (1982).
  * Kline and Zaffarano (1954) R. M. Kline and D. J. Zaffarano, Phys. Rev. 96, 1620 (1954).
  * Bienlein and Pleasonton (1962) J. Bienlein and F. Pleasonton, Nuclear Physics 37, 529 (1962).
  * Barker  _et al._ (1981) P. Barker, T. Ko, and M. Scandle, Nuclear Physics A 372, 45 (1981).


# Infinite Diameter Confidence Sets in Hedges’ Publication Bias Model

Jonas Moss
Department of Mathematics, University of Oslo
PB 1053, Blindern, NO-0316, Oslo, Norway __
<EMAIL_ADDRESS>
###### Abstract

Meta-analysis, the statistical analysis of results from separate studies, is a
fundamental building block of science. But the assumptions of classical meta-
analysis models are not satisfied whenever publication bias is present, which
causes inconsistent parameter estimates. Hedges’ selection function model
takes publication bias into account, but estimating and inferring with this
model is tough. Using a generalized Gleser–Hwang theorem, we show there is no
confidence set of guaranteed finite diameter for the parameters of Hedges’
selection model. This result provides a partial explanation for why inference
with Hedges’ selection model is fraught with difficulties.

_K_ eywords meta-analysis $\cdot$ confidence intervals $\cdot$ file-drawer
problem $\cdot$ publication bias $\cdot$ selection models $\cdot$ weight
function models

##  1 Introduction

A meta-analysis is a statistical analysis that quantitatively combines results
from separate scientific studies. When the studies measure the same
phenomenon, pooling of information allows us to predict the common effect size
with larger precision than we could have done with one study alone. Meta-
analyses are ubiquitous in the empirical sciences and forms a key component of
most systematic reviews such as Cochrane reviews (Higgins ., 2019).

Most meta-analytic techniques assume honest and unbiased reporting of results.
But there is ample evidence that the scientific literature is not unbiased, as
studies with significant p-values tend to be published with a greater
probability than other studies (Easterbrook ., 1991), a phenomenon called
publication bias by Sterling (1959) and the file-drawer problem by Rosenthal
(1979). When publication bias is present, there is no reason to trust the
results of meta-analytic methods that do not account for it, as the parameter
estimates will be inconsistent (Carter ., 2019).

Hedges’ (1992) publication bias model takes publication bias explicitly into
account using a selection model, and is arguably the most appropriate model
for publication bias (Carter ., 2019). Despite there being an $\mathtt{R}$ (R
Core Team, 2020) package for maximum likelihood estimation of this model,
called weightr (Coburn ., 2019), the model has not yet taken off. Its maximum
likelihood estimation methods are numerically unstable and its estimates can
be off even when they converge (Coburn ., 2019; Stanley, 2005). The estimate
of the mean effect size may be negative and of unrealistically large
magnitude, and the estimated heterogeneity parameter might be improbably
large. It turns out there are ridges in the likelihood function that can be
linked to this behavior (McShane ., 2016), but it has not been stated in clear
terms exactly what the consequences are for inferential procedures. The
purpose of this note is to explain why Hedges’ publication bias performs
poorly, by showing there is no confidence set for the mean effect size that
has infinite diameter with probability zero.

##  2 Hedges’ publication bias model

The most popular and well-known meta-analysis method is the random effects
model with normal likelihoods (Hedges Vevea, 1998). Written in hierarchical
notation, it equals

$\displaystyle\theta_{i}$ $\displaystyle\sim$ $\displaystyle
N(\theta_{0},\tau),$ $\displaystyle x_{i}\mid\theta_{i},\sigma$
$\displaystyle\sim$ $\displaystyle N(\theta_{i},\sigma_{i}).$

xzHere $x_{i}$ is the effect size and $\sigma_{i}$ is the standard deviation
of the $i$th study, $i=1,\ldots,N$. Following the convention in meta-analysis,
we assume all $\sigma_{i}$s to be known. The mean parameter $\theta_{0}$ is
the population effect size, $N(\theta_{0},\tau)$ is the effect size
distribution, and $\tau$ is the heterogeneity parameter. The purpose of the
effect size distribution is to model the fact that most effect size estimates
plugged into a meta-analysis do not appear to measure the same phenomenon. By
integrating out $\theta_{i}$, we find the density of $x_{i}$,

$f(x_{i};\theta_{0},\tau,\sigma_{i})=\phi(x_{i};\theta_{0},{(\tau^{2}+\sigma_{i}^{2})^{1/2}}),$

where $\phi$ is the density of a normal random variable.

We will assume that the random effects meta-analysis model is true in the
absence of publication bias. The mechanisms that cause publication bias modify
the density in a suitable way. Consider the case when only significant studies
at some specified level $\alpha$ are published. Assuming one-sided tests, the
p-values are $u_{i}=\Phi(-x_{i}/\sigma_{i})$, or normal one-sided p-values. We
will only deal with with one-sided p-values in this paper, as there is
usually, but not always, just one direction that is interesting to
researchers, reviewers, and editors. A one-sided p-value can also be used if
the researchers reported a two-sided value, since $p=0.05$ for a two-sided
hypothesis corresponds to $p=0.025$ for a one-sided hypothesis, et cetera.

Define $c_{\alpha}=\Phi^{-1}(1-\alpha)$, the cutoff for significance at level
$\alpha$. The _basic publication bias model_ is a truncated normal model with
density

$f(x_{i};\theta_{0},\sigma_{i})=\Phi\left(\frac{\theta_{0}-c_{\alpha}}{{(\sigma_{i}^{2}+\tau^{2})^{1/2}}}\right)^{-1}\phi(x;\theta_{0},{(\sigma_{i}^{2}+\tau^{2})^{1/2}})1(x_{i}/\sigma_{i}>c_{\alpha}],$
(2.1)

where $1[A]$ is the characteristic function of $A$. This model for publication
bias was introduced by Hedges (1984) in the context of $F$-distributions.

The basic publication bias model is unrealistic. It requires that no non-
significant studies are published. But even in the fields most severely
affected by publication bias, such as psychology, a non-negligible number of
non-significant studies are published (Motyl ., 2017). Moreover, the basic
publication bias model does not allow for different cutoffs for significance.
It is likely that some editors will accept studies reaching a significance at
$\alpha=0.025$, corresponding to $x_{i}/\sigma_{i}>1.96$ but not at
$\alpha=0.05$, corresponding to $x_{i}/\sigma_{i}>1.64$.

These problems can be rectified by adopting the selection model for
publication bias of Iyengar Greenhouse (1988), which models the following
scenario.

> Publication bias scenario. Alice the editor receives a study with the
> p-value $u$. Her publication decision is a random function of this p-value.
> That is, she will publish the study with some probability $w(u)$ and reject
> it with probability $1-w(u)$. Every study you will ever read in Alice’s
> journal has survived this selection mechanism, the rest are lost forever.

Let $w(u_{i})$ be a function of the p-value $u_{i}=\Phi(-x_{i}/\sigma_{i})$
taking values in $[0,1]$. Then $w(u_{i})$ is a probability for every $u_{i}$,
and the selection model

$f(x_{i};\theta_{0},{(\sigma_{i}^{2}+\tau^{2})^{1/2}})\propto\phi(x_{i};\theta_{0},{(\tau^{2}+\sigma_{i}^{2})^{1/2}})w(u)$
(2.2)

models the publication bias scenario exactly. This model can be viewed as a
rejection sampling procedure (Flury, 1990; von Neumann, 1951), where $\phi$
serves as proposal distribution for $f$. Variants of this model, with and
without covariates, has been studied by e.g. Dear Begg (1992); Vevea Hedges
(1995); Vevea Woods (2005); Citkowicz Vevea (2017).

Hedges (1992) studies the selection model when $w$ is a step function with
fixed steps. Let $\alpha$ be a vector with elements
$0=\alpha_{0}<\alpha_{1}<\ldots<\alpha_{K}=1$ and $\rho$ be a $K$-ary non-
negative, non-decreasing vector having its first element equal to $\rho_{1}=1$
for identifiability. Define the step function $w$ based on $\alpha$ and $\rho$
as

$w(u;\rho,\alpha)=\sum_{k=1}^{K}\rho^{k}1_{(\alpha_{k-1},\alpha_{k}]}(u).$
(2.3)

We call the selection model with a step function _Hedges’ publication bias
model_. Its density is

$f(x_{i};\theta_{0},{(\tau^{2}+\sigma_{i}^{2})^{1/2}})\propto\sum_{k=1}^{K}\rho^{k}\phi(x_{i};\theta_{0},{(\tau^{2}+\sigma_{i}^{2})^{1/2}})1{(\alpha_{k-1},\alpha_{k}]}(u_{i}).$
(2.4)

Interpreting Hedges’ publication bias model is easy. When the editor receives
a study with p-value $u$, she finds the $k$ such that
$u\in(\alpha_{k-1},\alpha_{k}]$ and accepts with probability $\rho^{k}$. Since
$\rho^{1}=1$, she always accepts when $u\in[0,\alpha_{1}]$. The vector $\rho$
is non-increasing since a publication decision based solely on p-values should
always act favorably towards lower p-values. The parameters
$(\theta_{0},\tau,\rho)$ of the model are identified when $\alpha$ is fixed
(Moss De Bin, 2020, Web Appendix A).

Hedges’ publication bias model allows both for non-significant studies to be
published and allows the editor to act differently towards different cutoffs
such as $\alpha=0.025$ and $\alpha=0.05$. In addition, the model can
approximate any non-increasing selection function $w$ by increasing the number
of steps.

We can write Hedges’ model as a mixture model on the form

$f(x_{i};\theta_{0},{(\tau^{2}+\sigma_{i}^{2})^{1/2}})=\sum_{k=1}^{K}\pi^{k}f^{k}(x_{i};\theta_{0},{(\tau^{2}+\sigma_{i}^{2})^{1/2}}).$
(2.5)

Where $f^{k},k\leq K$ are normal densities truncated to
$(\Phi^{-1}(1-\alpha_{k-1}),\Phi^{-1}(1-\alpha_{k})]$, and $\pi^{k}$ are
mixture probabilities, i.e., $\pi^{k}>0$ for each $k$ and
$\sum_{k=1}^{K}\pi^{k}=1$. The mixture probabilities $\pi^{k}$ are functions
of $(\theta_{0},\tau,\sigma_{i},\rho)$, see the appendix (p. Appendix) for
their formula.

The main benefit of Hedges’ publication bias model (2.4) is how it models
p-values based publication bias directly, there is no approximation involved.
If you believe in the random effects meta-analysis model and the p-value based
publication bias scenario, Hedges’ publication bias model is simply the
correct model. Most statistical methods correcting for publication bias in the
literature either do not make use of an explicit statistical model or do not
estimate the parameters $\theta_{0}$ and $\tau$. For instance, the funnel plot
of Egger . (1997) is a graphical method, while the trim-and-fill method of
Duval Tweedie (2000) is a non-parametric method based on making the funnel
plot symmetric. Stanley (2005); Stanley Doucouliagos (2014) discuss various
misspecified regression-based estimators of the corrected effect size $\theta$
based on the fixed effect Hedges’ publication bias model. The estimating
p-curve method of Simonsohn . (2014) and p-uniform of van Assen . (2015); van
Aert . (2016) are two methods for dealing with publication bias hailing from
psychology. Both are based on a variant of the basic publication bias model,
but with fixed instead of random effects, and both employ somewhat unusual
estimation methods (McShane ., 2016). Since there is ample evidence of
heterogeneity in meta-analysis, restricting oneself to the fixed effects meta-
analysis is a mistake.

Hedges’ model has some downsides. It models only bias due to selection of
p-values, not every source of bias, such as language bias (Egger Smith, 1998).
Second, it may not be best model for biases with other causes than the
publication process, such as p-hacking (Simmons ., 2011). Moss De Bin (2020)
propose a related model that may be more successful at correcting for
p-hacking.

Estimation of Hedges’ model is difficult. Stanley (2005, Section 6.3)
discusses three problematic cases in economics when maximum likelihood was
used to estimate Hedges’ publication bias model. McShane . (2016, Appendix A)
notes that while estimation of the basic publication bias model is hard,
introducing the heterogeneity parameter exacerbates the problem. The
likelihood function has contours following approximately
$\tau\propto|\theta|^{1/2}$.

Figure 2.1: Contour lines for the log-likelihood for the simple publication
model using power posing data of Cuddy . (2018).

Figure 2.1 shows the contour lines for the meta-analysis of Cuddy . (2018)
where the selection probabilities of the step function (2.3) are fixed at
$\rho=(1,0.6,0.1)$ and $\alpha=(0,0.025,0.05,1$. Only just around the maximum
at $\hat{\theta}_{0}=0.55$ and $\hat{\sigma}=5\cdot 10^{-7}$ can the
likelihood be approximated with a quadratic function.

##  3 Confidence sets of infinite diameter

Fix some measurable space $(\Omega,\mathcal{F})$, let $\mathcal{P}$ be a
family of dominated probability measures defined on this measurable space, and
$\Pi$ a partition of $\mathcal{P}$. Recall that a partition of $\mathcal{P}$
is a collection of disjoint non-empty subsets $\pi$ of $\mathcal{P}$ such that
$\bigcup_{\pi\in\Pi}\pi=\mathcal{P}$. When $p$ is a density associated with a
$P\in\mathcal{P}$, we will use the standard notation $[p]$ to denote the
unique part $\pi$ containing $P$.

######  Definition 1.

A _confidence set_ of level $\alpha$ is a family of rejection sets
$\\{R(\pi)\\},\pi\in\Pi$ such that

$\sup_{\pi\in\Pi}\sup_{P\in\pi}P(R(\pi))\leq\alpha.$

If the inequality is an equality, the confidence set has _size_ $\alpha$.

This definition of confidence sets might look slightly unfamiliar, but it is a
straight-forward generalization of the definitions in Casella Berger (2002,
Definition 9.1.5) and Lehmann Romano (2005, Section 3.5). Usually, a
confidence is defined as a set $C$ adhering to the relation

$\pi_{0}\in C\iff\omega\notin R(\pi_{0}).$ (3.1)

That is, $\pi_{0}\in C$ if and only if we accept the null-hypothesis
$H_{0}:\pi=\pi_{0}$, but we will only need the formulation using rejection
sets in this paper. When confidence sets are defined in terms of rejection
sets, there is sometimes no partition $\Pi$ to take the supremum over, and the
definition reduces to $\sup_{P\in\mathcal{P}}P(R(P))\leq\alpha$. The term
_confidence interval_ is far more common than confidence set, but this
requires that the $C$ in equation 3.1 is an interval, which we will not
require here.

The following example should make Definition 1 familiar.

######  Example 2.

Consider the usual _t_ -confidence interval with $n$ observations. In this
case, $\mathcal{P}$ contains all measures $P_{\mu,\sigma}^{n}$, where
$P_{\mu,\sigma}$ is the probability measure of a normal with mean $\mu$ and
standard deviation $\sigma$. Here $Q^{n}$ denotes the $n$-fold product measure
of $Q$, corresponding to $n$ independent samples from $Q$ when $Q$ is a
probability measure. Let $\pi(\mu)=\\{P_{\mu,\sigma}^{n}\mid\sigma>0\\}$
contain all normal probability measures with mean $\mu$ and some positive
standard deviation. Then $\\{\pi(\mu)\\},\mu\in\mathbb{R}$ defines a partition
of $\mathcal{P}$. A two-sided _t_ -confidence interval is a confidence set of
size $\alpha$ with partition $\\{\pi(\mu)\\},\mu\in\mathbb{R}$ according to
Definition 1.

Now we must find out how to measure the size of confidence sets. Let a _size
function_ be any positive function $\left\|\cdot\right\|:\Pi\to[0,\infty)$. In
the _t_ -confidence interval example above, $\left\|\cdot\right\|$ can be
taken to be $\left\|\pi\right\|=|\mu|$ for the unique $\mu$ associated with
each $\pi$.

The diameter of a confidence set is the random variable

$D(\omega)=\sup_{\pi\in\Pi}\\{\left\|\pi\right\|\mid\omega\notin R(\pi)\\}.$
(3.2)

The diameter tells you the size of the largest accepted $\pi$. We will assume
that $D$ is Borel measurable.

######  Definition 3.

A confidence set has infinite diameter with $P$-positive probability if
$P(D=\infty)>0$. It has infinite diameter with positive probability if
$P(D=\infty)>0$ for all $P\in\mathcal{P}$.

This is a modified Gleser–Hwang theorem (Gleser Hwang, 1987, Theorem 1)
formulated in terms of partitions. Its proof is in the appendix (p. Appendix).

######  Theorem 4 (Gleser–Hwang theorem).

Suppose there is a sequence $\\{p_{n}\\}$ of densities derived from
$\mathcal{P}$ satisfying the following:

  1. (i)

There is a density $p^{\star}$ such that $p_{n}$ converges to $p^{\star}$
pointwise,

  2. (ii)

$\operatorname{supp}p\supseteq\operatorname{supp}p^{\star}$ for all densities
$p$ derived from $\mathcal{P}$,

  3. (iii)

the size of the equivalence class $[p_{n}]$ goes to infinity as $n$ increases,
$\left\|[p_{n}]\right\|\to\infty$.

Then every confidence set with level $\alpha>0$ has infinite diameter with
positive probability.

Following the terminology of Berger . (1999), we will say that families
$\mathcal{P}$ of probabilities satisfying the conclusion of Theorem 4 for a
suitable partition $\Pi$ belong to the Gleser–Hwang class. To make the
Gleser–Hwang class more familiar, we will present two examples. More examples
can be found in the papers of Gleser Hwang (1987) and Berger . (1999).

Fieller’s problem is the best known case of a badly behaved confidence set.
Let $(X,Y)$ be an observation from a bivariate normal
$N([\mu_{1},\mu_{2}],I\sigma^{2})$, where $\sigma^{2}$ is known. We want to
form a confidence set for the ratio $\mu_{2}/\mu_{1}$. The most famous
confidence set is due to Fieller (1940). His confidence set can be finite, the
whole real line, or the union of two disjoint semi-infinite intervals, all
with positive probability (Koschat, 1987).

But it is not only Fieller’s confidence set that might be infinitely long. The
Gleser–Hwang theorem can be used to show that every confidence set for
$\mu_{2}/\mu_{1}$ must be infinitely long with positive probability. This
result is almost independent of the distribution of $X$ and $Y$. To state this
result in our notation, let $\mathcal{P}$ be a family of bivariate
distributions over $(X,Y)$. All of these distributions have the same support,
and all of them have finite means $E_{p}(X)$ and $E_{p}(Y)$. Moreover, assume
$E_{p}(X)=0$ is attainable for some $p\in\mathcal{P}$. Define the partition
$\Pi$ by $p,q\in\pi$ if and only if $E_{p}(X)/E_{p}(Y)=E_{q}(X)/E_{q}(Y)$, and
let $\left\|[p]\right\|=|E_{q}(X)/E_{q}(Y)|$, i.e., the ratio of means. Choose
a sequence $p_{n}(x,y)=p(x,y)$, where $p(x,y)$ is density with means
$E_{p}X>0$ and $E_{p}Y=0$. Then $\left\|[p]\right\|=\infty$, the conditions of
Theorem 4 are satisfied, and every confidence set with level $\alpha>0$ has
infinite diameter with positive probability.

Another example is due to Bahadur Savage (1956), who studies non-parametric
testing of the mean, and concludes the mean cannot be meaningfully tested.
They are working with a family $\mathcal{P}$ of densities over $\mathbb{R}$
that covers all finite means, has finite variances, and is closed under convex
combinations. Similar problems were considered by Romano (2004) and Donoho
(1988).

Using the Gleswer–Hwang theorem, it is easy to verify that every confidence
set has infinite diameter with positive probability. Define the partition
$\Pi$ by $p,q\in\pi$ if and only if $E_{p}(X)=E_{q}(X)$, and let
$\left\|[p]\right\|=|E_{p}(X)|$. Let

$p_{n}(x)=\left(1-\frac{1}{n}\right)q_{0}(x)+\frac{1}{n}q_{n^{2}}(x),$

where $q_{0}$ has mean $0$ and $q_{n^{2}}$ has mean $n^{2}$. Then
$\left\|[p_{n}]\right\|=n$, the conditions of Theorem 4 are satisfied, and
every confidence set with level $\alpha>0$ has infinite diameter with positive
probability.

There are several natural candidates for $\Pi$ when working with Hedges’
selection function model (2.4). We will work with three of them. First,
consider the partition where $p,q\in\pi$ if and only if $p,q$ have the same
mean effect size parameter $\theta_{0}$. We will equip this partition with the
size function $\left\|\pi\right\|=|\theta_{0}|$, and it corresponds to a
confidence set for $\theta_{0}$. Second, consider the partition where all
$p,q\in\pi$ have the same heterogeneity parameter $\tau$, equipped with
$\left\|\pi\right\|=\tau$. Finally, we will work with the partition where all
$p,q\in\pi$ have the same heterogeneity parameter $\tau$ and population effect
size $\theta_{0}$, and equip it with
$\left\|\pi\right\|={(\theta_{0}^{2}+\tau^{2})^{1/2}}$. This information is
summarized in Table 1 for convenience.

Table 1: The three partitions $\Pi$ for the selection function model | Symbol | Size $\left\|\cdot\right\|$ | Confidence set
---|---|---|---
Mean effect size | $\theta_{0}$ | $\left\|\pi\right\|=|\theta_{0}|$ | Confidence set for $\theta_{0}$
Heterogeneity parameter | $\tau$ | $\left\|\pi\right\|=\tau$ | Confidence set for $\tau$
Both parameters | $(\theta_{0},\tau)$ | $\left\|\pi\right\|={(\theta_{0}^{2}+\tau^{2})^{1/2}}$ | Joint confidence set for $(\theta,\tau)$

Let us take a look at the basic publication bias model (2.1) again. To use
Theorem 4 we need a witnessing sequence of functions $p_{n}\to p$ satisfying
the conditions (ii) and (iii). The next lemma shows how to make such a witness
for the truncated normal. Its proof is in the appendix (p. Appendix).

######  Lemma 5.

Let $f_{n}$ be a normal density truncated to $[a,b)$, where $b=\infty$ is
allowed, with underlying mean $\theta_{n}=-n$ and standard deviation
$\sigma_{n}^{2}=n+c$ for some $c\in\mathbb{R}$. Then $f_{n}$ converges
pointwise to $\exp(-x)/[\exp(-a)-\exp(-b)]$, the distribution of an
exponential variable truncated to $[a,b)$.

Using Lemma 5 it is not hard to show that the basic publication bias model
(2.1) is a member of the Gleser–Hwang class.

######  Theorem 6.

Assume we have $N$ independent samples from the basic publication bias model
(2.1). Then any confidence set for $\theta_{0},\tau$, or $(\theta_{0},\tau)$
with level $\alpha>0$ will have infinite diameter with positive probability.

###### Proof.

Let $\Pi$ be the partition of $\mathcal{P}$ where $p,q\in\pi$ if and only if
they share the same $\theta_{0}$, and let $||[p]||=|\theta_{0}|$. We are
dealing with products of densities of the form (2.1), that is,

$p(x)=\prod_{i=1}^{N}\Phi\left(\frac{\theta_{0}-c_{\alpha}}{{(\sigma_{i}+\tau)^{1/2}}}\right)^{-1}\phi(x_{i};\theta_{0},{(\sigma_{i}^{2}+\tau^{2})^{1/2}}),$

where $\sigma_{i}$ are known parameters. From Lemma 5, $p_{n}$ converges to a
product of truncated exponentials when $\theta_{n}=-n$ and $\tau_{n}^{2}=n$.
Since $||[p_{n}]||=n$, the three conditions of Theorem 4 are satisfied. The
proofs for $||[p]||=\tau$ and $||[p]||={(\tau^{2}+\theta_{0}^{2})^{1/2}}$ are
similar and omitted. ∎

Proving the analogue of Theorem 6 for Hedges’ publication is only somewhat
more involved. We will use the mixture representation of (2.5) and a lemma
generalizing Theorem 4 to a certain kind of mixtures.

Let $f^{1},f^{2},\ldots,f^{K}$ be a sequence of densities,
$\pi=(\pi^{1},\pi^{2},\ldots,\pi^{K})$ be a probability vector, and
$p=\sum_{k\leq K}\pi^{k}f^{k}$ be a mixture distribution. We will assume that
the size of $[p]$ equals the size of any of its mixture components $[f^{k}]$
for some size $\left\|\cdot\right\|$, i.e.,
$\left\|[p]\right\|=\left\|[f^{k}]\right\|$ for all $k$. Why we do this will
be clear in the proof of Theorem 8, but think of it this way: If all of $p$s
mixture components have the same mean, the mean of $p$ equals the mean of any
$f^{k}$.

######  Lemma 7.

Let $\mathcal{P}$ be a class of $K$-ary mixture distributions and $||[p]||$ be
as assumed above. Assume there is a sequence $p_{n}=\sum_{k\leq
K}\pi_{n}^{k}f_{n}^{k}$ and a subset $K^{\prime}$ such that

  1. (i)

For all $k\in K^{\prime}$, there is a density $f^{k\star}$ such that
$f_{n}^{k}$ converges to $f^{k\star}$ pointwise.

  2. (ii)

For all mixtures $p$,
$\operatorname{supp}p\supseteq\operatorname{supp}f^{k\star}$ for all $k\in
K^{\prime}$.

  3. (iii)

For all $k\in K^{\prime}$, the size of $[f_{n}^{k}]$ goes to infinity,
$\left\|[f_{n}^{k}]\right\|\to\infty$.

  4. (iv)

The density concentrates on the components indexed by $K^{\prime}$,
$\lim_{n\to\infty}\sum_{k\in K^{\prime}}\pi_{n}^{k}=1$.

Then every confidence set with level $\alpha>0$ has infinite diameter with
positive probability.

###### Proof.

We employ Theorem 4. By (i) and (iv), $p_{n}$ converges pointwise to the
density

$\sum_{k\in K^{\prime}}\pi_{n}^{k\star}f_{n}^{k\star}=p^{\star}.$

That $\operatorname{supp}p\supseteq p^{\star}$ follows from (ii) and (iv).
Finally, from the assumption
that$\left\|[p_{n}]\right\|=\left\|[f_{n}^{k}]\right\|$, we get that
$\left\|[p_{n}]\right\|\to\infty$ too. ∎

######  Theorem 8.

Assume we have $N$ independent samples from the publication bias model (2.1),
where the selection probability $\rho$ is unknown and $\alpha$ is known. Then
any confidence set for $\theta_{0},\tau$, or $(\theta_{0},\tau)$ will have
infinite diameter with positive probability.

###### Proof.

Let $n=1$ and consider confidence sets for $\theta_{0}$. Let $\Pi$ be the
partition of $\mathcal{P}$ where $p,q\in\pi$ if and only if they share the
same $\theta_{0}$, and let $||[p]||=|\theta_{0}|$. Then
$||f^{k}||=|\theta_{0}|$ from the mixture representation (2.5). Using Lemma 5,
we see that $f^{k},k<K$ converges pointwise to truncated exponentials when
$\theta_{0}=-n$ and $\tau_{n}^{2}=n$, so that (i), (iii) of Proposition 7 are
satisfied with the set $K^{\prime}=\\{1,2,\ldots,K-1\\}$. Moreover, since we
assume that $\rho$ is decreasing, (ii) is satisfied as well. The mixture
probabilities for $k\neq K$ can be fixed at e.g. $\pi=1/(K-1)$, and (iv) is
satisfied as well. The proofs for $||[p]||=\tau$ and
$||[p]||={(\tau^{2}+\theta_{0}^{2})^{1/2}}$ are similar and omitted.

When $N>1$, expand the expression
$\prod_{i=1}^{N}\sum_{k<K}\pi_{i}^{k}(\sigma_{i})f_{i}^{k}(\sigma_{i}),$ and
use the same reasoning as in the first part of this proof. ∎

##  4 Remarks

Well-behaved confidence sets for Hedges publication bias model do not exist,
but well-behaved credibility sets do. Bayesian estimation of Hedges’ model can
be made routine, as it is easy to find uncontroversial priors for $\theta_{0}$
and $\tau$. In practical meta-analyses we know that $\theta_{0}$ cannot be
large, and is likely to be close to $0$. Moreover, since it is common effects
in meta-analyses to be interpreted as the aggregation of many small effects,
the central limit theorem justifies using a normal prior. As we want to remove
prior mass from negative $\theta_{0}$s of large magnitude, $N(0,1)$ is a
decent standard prior. Similarly, a half-normal is a reasonable prior for the
heterogeneity parameter $\tau$. Moss De Bin (2020) employed these priors on
several examples.

## Appendix

The following sandwich convergence theorem is used in the proof of the
Gleser–Hwang theorem.

######  Lemma 9 (Billingsley (1995, Exercise 16.4(a))).

Suppose $a_{n},b_{n},f_{n}$ converge to $a,b,f$ and $a_{n}\leq f_{n}\leq
b_{n}$ for all $n$. If $\int a_{n}d\mu\to\int ad\mu$ and $\int
b_{n}d\mu\to\int bd\mu$, then $\int f_{n}d\mu\to\int fd\mu$ for any measure
$\mu$.

The proof of Theorem 4 closely follows the proof of Gleser Hwang (1987,
Theorem 1).

###### Proof of Theorem 4.

We can assume without loss of generality that $\left\|[p_{n}]\right\|\geq n$,
as we can choose a suitable sub-sequence if we have to. By definition of the
diameter $D$ (3.2) we see that

$\\{D\geq n\\}=\\{\omega\in\Omega\mid\textrm{there is a }\pi\textrm{ such
that}\left\|\pi\right\|\geq n\textrm{ and }\omega\in R^{c}(\pi)\\}.$

It follows that, if $||[p_{n}]||\geq n$, then $R^{c}([p_{n}])\subseteq\\{D\geq
n\\}$. Since we assume that $\left\|[p_{n}]\right\|\geq n$ and

$1-\alpha\leq P_{n}(R^{c}([p_{n}]))=\int_{R^{c}([p_{n}])}p_{n}d\mu$

by definition of a confidence set, we have that

$0<1-\alpha\leq\int_{R^{c}([p_{n}])}p_{n}d\mu\leq\int_{D\geq n}p_{n}d\mu$
(4.1)

for all $n$. Since $p_{n}$ and $p^{\star}$ are densities,

$\lim_{n\to\infty}\int p_{n}d\mu=1=\int
p^{\star}d\mu=\int\lim_{n\to\infty}p_{n}d\mu.$

This allows us to use Lemma 9 with $a_{n}=0$, $b_{n}=p_{n}$, and
$f_{n}=1_{D\geq n}p_{n}$ and conclude that

$\int_{D\geq n}p_{n}d\mu\to\int_{D=\infty}p^{\star}d\mu.$ (4.2)

Combining equations (4.1) and (4.2), we get

$0<1-\alpha\leq\int_{D=\infty}p^{\star}d\mu.$

Let $P\in\mathcal{P}$ be arbitrary, $p$ be its density, and consider

$P(D=\infty)=\int_{D=\infty}pd\mu\geq\int_{D=\infty\cap\operatorname{supp}p^{\star}}\left(\frac{p}{p^{\star}}\right)p^{\star}d\mu.$

Since $\int_{D=\infty}p^{\star}d\mu>0$ and $p/p^{\star}>0$ on
$\operatorname{supp}p^{\star}$ (since
$\operatorname{supp}p\supseteq\operatorname{supp}p^{\star}$ by assumption), we
see that
$\int_{D=\infty\cap\operatorname{supp}p^{\star}}(p/p^{\star})p^{\star}d\mu>0$
too. It follows that $P(D=\infty)>0$, and, since $P$ is arbitrary, $D$ has
infinite diameter with positive probability. ∎

Now we prove Lemma 5.

###### Proof of Lemma (5).

Let $n>-c$, so that $\sigma_{n}^{2}>0$. Recall the well-known formula for the
normal truncated to $[a,b]$, and substitute $\theta_{n}=-n$ and
$\sigma^{2}_{n}=n+c$,

$\displaystyle f_{n}(x)$ $\displaystyle=$
$\displaystyle\frac{1}{\Phi\left(\frac{b-\theta_{n}}{\sigma_{n}}\right)-\Phi\left(\frac{a-\theta_{n}}{\sigma_{n}}\right)}\phi(x;\theta_{n},\sigma_{n})1[a,b](x),$
$\displaystyle=$
$\displaystyle\frac{\phi(x;-n,(n+c)^{1/2})1[a,b](x)}{\Phi[-(a+n)(n+c)^{-1/2}]-\Phi[-(b+n)(n+c)^{-1/2}]}.$

The normal density part equals

$\phi(x;-n,(n+c)^{1/2})=(2\pi)^{-1/2}(n+c)^{-1/2}\exp\left(-\frac{x^{2}+2nx+n^{2}}{2(n+c)}\right).$

When $n$ is large compared to $x$, the term $x^{2}/2(n+c)$ is negligible,
hence

$\displaystyle\phi(x;-n,(n+c)^{1/2})$ $\displaystyle\approx$
$\displaystyle(2\pi)^{-1/2}(n+c)^{-1/2}\exp(-n^{2}/2(n+c))\exp(-x),$
$\displaystyle=$ $\displaystyle(n+c)^{-1/2}\phi(n/(n+c)^{1/2})\exp(-x).$

From Equation 5 of Borjesson Sundberg (1979), we know that
$\Phi(-x)\approx\phi(x)/x$ as $x$ grows. Then

$\Phi[-(a+n)(n+c)^{-1/2}]\approx\frac{(n+c)^{1/2}}{n+a}\phi[-(a+n)(n+c)^{-1/2}],$

and using the same reasoning as above, we find that
$\phi[-(a+n)(n+c)^{-1/2}]\approx\phi((n+c)^{1/2})\exp(-a)$ as $n$ increase.
Therefore,

$\Phi[-(a+n)(n+c)^{-1/2}]\to\frac{(n+c)^{1/2}\phi((n+c)^{1/2})\exp(-a)}{a+n}.$

Since this reasoning applies to $b$ as well, we get that $f$ approaches

$\displaystyle\frac{(n+c)^{-1/2}\phi(n/(n+c)^{1/2})\exp(-x)}{\Phi[-(a+n)n^{-1/2}]-\Phi[-(b+n)n^{-1/2}]},$
$\displaystyle\approx$
$\displaystyle\frac{(n+c)^{-1/2}\phi(n/(n+c)^{1/2})\exp(-x)}{(n+c)^{1/2}\phi((n+c)^{1/2})\left[\frac{\exp(-a)}{a+n}-\frac{\exp(-b)}{b+n}\right]},$
$\displaystyle=$
$\displaystyle\frac{\phi(-n^{2}/2(n+c))\exp(-x)}{\phi((n+c)^{1/2})\left[\frac{\exp(-a)}{a+n}-\frac{\exp(-b)}{b+n}\right]},$
$\displaystyle\approx$
$\displaystyle\exp(-x)n^{-1}\left[\frac{\exp(-a)}{a+n}-\frac{\exp(-b)}{b+n}\right]^{-1},$
$\displaystyle\to$ $\displaystyle\exp(-x)\left[\exp(-a)-\exp(-b)\right]^{-1}.$

Here the third line follows from $\phi(n/(n+c)^{1/2})/\phi((n+c)^{1/2})\to 1$.
∎

These are the formula for the mixture probabilities $\pi_{i}$, see (2.5). Let
$c_{k}=\Phi^{-1}(1-\alpha_{k})$ and define

$c=\sum_{k=1}^{K}\rho^{k}[\Phi(c_{k-1};\theta_{0},(\tau^{2}+\sigma_{i}^{2})^{1/2})-\Phi(c_{k};\theta_{0},(\tau^{2}+\sigma_{i}^{2})^{1/2})].$

Then

$\pi^{k}=c^{-1}\rho^{k}[\Phi(c_{k-1};\theta_{0},(\tau^{2}+\sigma_{i}^{2})^{1/2})-\Phi(c_{k};\theta_{0},(\tau^{2}+\sigma_{i}^{2})^{1/2})].$

## References

  * Bahadur Savage (1956) bahadur1956nonexistenceBahadur, RR. Savage, LJ.  1956\.  The Nonexistence of Certain Statistical Procedures in Nonparametric Problems The nonexistence of certain statistical procedures in nonparametric problems.  Annals of Mathematical Statistics2741115–1122.  https://doi.org/10.1214/aoms/1177728077
  * Berger . (1999) berger1999integratedBerger, JO., Liseo, B. Wolpert, RL.  1999\.  Integrated Likelihood Methods for Eliminating Nuisance Parameters Integrated likelihood methods for eliminating nuisance parameters.  Statistical Science1411–28.  https://doi.org/10.1214/ss/1009211804
  * Billingsley (1995) billingsley1995probabilityBillingsley, P.  1995\.  Probability and Measure Probability and measure.  John Wiley & Sons.
  * Borjesson Sundberg (1979) borjesson1979simpleBorjesson, P. Sundberg, C.  1979\.  Simple Approximations of the Error Function Q(x) for Communications Applications Simple approximations of the error function Q(x) for communications applications.  IEEE Transactions on Communications273639–643.  https://doi.org/10.1109/TCOM.1979.1094433
  * Carter . (2019) Carter2019-rwCarter, EC., Schönbrodt, FD., Gervais, WM. Hilgard, J.  2019\.  Correcting for Bias in Psychology: A Comparison of Meta-Analytic Methods Correcting for bias in psychology: A comparison of meta-analytic methods.  Advances in Methods and Practices in Psychological Science22115–144.  https://doi.org/10.1177/2515245919847196
  * Casella Berger (2002) Casella2002-lgCasella, G. Berger, RL.  2002\.  Statistical Inference Statistical inference.  Duxbury Press.
  * Citkowicz Vevea (2017) citkowicz2017parsimoniousCitkowicz, M. Vevea, JL.  2017\.  A Parsimonious Weight Function for Modeling Publication Bias A parsimonious weight function for modeling publication bias.  Psychological Methods22128–41.  https://doi.org/10.1037/met0000119
  * Coburn . (2019) Coburn2019-ecCoburn, KM., Vevea, JL. Coburn.  2019\.  R package ‘weightr’ R package ‘weightr’ [].  cran.rstudio.org.  https://cran.rstudio.org/web/packages/weightr/weightr.pdf
  * Cuddy . (2018) Cuddy2018-kpCuddy, AJC., Schultz, SJ. Fosse, NE.  2018\.  P-Curving a More Comprehensive Body of Research on Postural Feedback Reveals Clear Evidential Value for Power-Posing Effects: Reply to Simmons and Simonsohn (2017) P-Curving a more comprehensive body of research on postural feedback reveals clear evidential value for power-posing effects: Reply to Simmons and Simonsohn (2017).  Psychological Science294656–666.  https://doi.org/10.1177/0956797617746749
  * Dear Begg (1992) Dear1992-gwDear, KBG. Begg, CB.  1992\.  An Approach for Assessing Publication Bias Prior to Performing a Meta-analysis An approach for assessing publication bias prior to performing a Meta-analysis.  Statistical Science72237–245.  https://doi.org/10.1214/ss/1177011363
  * Donoho (1988) Donoho1988-hgDonoho, DL.  1988\.  One-Sided Inference about Functionals of a Density One-sided inference about functionals of a density.  Annals of Statistics1641390–1420.  https://doi.org/10.1214/aos/1176351045
  * Duval Tweedie (2000) Duval2000-ctDuval, S. Tweedie, R.  2000\.  A Nonparametric “Trim and Fill” Method of Accounting for Publication Bias in Meta-analysis A nonparametric “trim and fill” method of accounting for publication bias in meta-analysis.  Journal of the American Statistical Association9544989–98.  https://doi.org/10.1080/01621459.2000.10473905
  * Easterbrook . (1991) Easterbrook1991-phEasterbrook, PJ., Berlin, JA., Gopalan, R. Matthews, DR.  1991\.  Publication Bias in Clinical Research Publication bias in clinical research.  The Lancet3378746867–872.  https://doi.org/10.1016/0140-6736(91)90201-y
  * Egger . (1997) Egger1997-rdEgger, M., Davey Smith, G., Schneider, M. Minder, C.  1997\.  Bias in Meta-analysis Detected by a Simple, Graphical Test Bias in meta-analysis detected by a simple, graphical test.  British Medical Journal3157109629–634.  https://doi.org/10.1136/bmj.315.7109.629
  * Egger Smith (1998) Egger1998-kjEgger, M. Smith, GD.  1998\.  Bias in Location and Selection of Studies Bias in location and selection of studies.  British Medical Journal316712461–66.  https://doi.org/10.1136/bmj.316.7124.61
  * Fieller (1940) Fieller1940-lgFieller, EC.  1940\.  The Biological Standardization of Insulin The biological standardization of insulin.  Supplement to the Journal of the Royal Statistical Society711–64.  https://doi.org/10.2307/2983630
  * Flury (1990) flury1990acceptanceFlury, BD.  1990\.  Acceptance–Rejection Sampling Made Easy Acceptance–rejection sampling made easy.  SIAM Review323474–476.  https://doi.org/10.1137/1032082
  * Gleser Hwang (1987) gleser1987nonexistenceGleser, LJ. Hwang, JT.  1987\.  The Nonexistence of 100(1 - $\alpha$)% Confidence Sets of Finite Expected Diameter in Errors-in-Variables and Related Models The nonexistence of 100(1 - $\alpha$)% confidence sets of finite expected diameter in errors-in-variables and related models.  Annals of Statistics1541351–1362.  https://doi.org/10.1214/aos/1176350597
  * Hedges (1984) hedges1984estimationHedges, LV.  1984\.  Estimation of Effect Size under Nonrandom Sampling: The Effects of Censoring Studies Yielding Statistically Insignificant Mean Differences Estimation of effect size under nonrandom sampling: The effects of censoring studies yielding statistically insignificant mean differences.  Journal of Educational and Behavioral Statistics9161–85.  https://doi.org/10.3102/10769986009001061
  * Hedges (1992) hedges1992modelingHedges, LV.  1992\.  Modeling Publication Selection Effects in Meta-analysis Modeling publication selection effects in meta-analysis.  Statistical Science72246–255.  https://doi.org/10.1214/ss/1177011364
  * Hedges Vevea (1998) hedges1998fixedHedges, LV. Vevea, JL.  1998\.  Fixed- and Random-Effects Models in Meta-analysis Fixed- and random-effects models in meta-analysis.  Psychological Methods34486.  https://doi.org/10.1037/1082-989x.3.4.486
  * Higgins . (2019) Higgins2019-vvHiggins, JPT., Thomas, J., Chandler, J., Cumpston, M., Li, T., Page, MJ. Welch, VA.  2019\.  Cochrane Handbook for Systematic Reviews of Interventions Cochrane handbook for systematic reviews of interventions.  John Wiley & Sons.
  * Iyengar Greenhouse (1988) iyengar1988selectionIyengar, S. Greenhouse, JB.  1988\.  Selection Models and the File Drawer Problem Selection models and the file drawer problem.  Statistical Science31109–117.  https://doi.org/10.1214/ss/1177013012
  * Koschat (1987) Koschat1987-dkKoschat, MA.  1987\.  A Characterization of the Fieller Solution A characterization of the fieller solution.  Annals of Statistics151462–468.  https://doi.org/10.1214/aos/1176350282
  * Lehmann Romano (2005) lehmann2006testingLehmann, EL. Romano, JP.  2005\.  Testing Statistical Hypotheses Testing statistical hypotheses.  Springer Science & Business Media.  https://doi.org/10.1007/0-387-27605-x
  * McShane . (2016) McShane2016-rbMcShane, BB., Böckenholt, U. Hansen, KT.  2016\.  Adjusting for Publication Bias in Meta-analysis: An Evaluation of Selection Methods and Some Cautionary Notes Adjusting for publication bias in Meta-analysis: An evaluation of selection methods and some cautionary notes.  Perspectives on Psychological Science115730–749.  https://doi.org/10.1177/1745691616662243
  * Moss De Bin (2020) moss2019modellingMoss, J. De Bin, R.  2020\.  Modelling Publication Bias and p-hacking Modelling publication bias and p-hacking.  arXiv preprints.  https://arxiv.org/abs/1911.12445
  * Motyl . (2017) Motyl2017-dxMotyl, M., Demos, AP., Carsel, TS., Hanson, BE., Melton, ZJ., Mueller, AB.Skitka, LJ.  2017\.  The State of Social and Personality Science: Rotten to the Core, Not So Bad, Getting Better, or Getting Worse? The state of social and personality science: Rotten to the core, not so bad, getting better, or getting worse?  Journal of Personality and Social Psychology113134–58.  https://doi.org/10.1037/pspa0000084
  * R Core Team (2020) Team2013-ttR Core Team.  2020\.  R: A Language and Environment for Statistical Computing R: A language and environment for statistical computing [].  Vienna, Austria.
  * Romano (2004) romano2004nonRomano, JP.  2004\.  On Non-parametric Testing, the Uniform Behaviour of the t-test, and Related Problems On non-parametric testing, the uniform behaviour of the t-test, and related problems.  Scandinavian Journal of Statistics314567–584.  https://doi.org/10.1111/j.1467-9469.2004.00407.x
  * Rosenthal (1979) Rosenthal1979-pmRosenthal, R.  1979\.  The File Drawer Problem and Tolerance for Null Results The file drawer problem and tolerance for null results.  Psychological Bulletin863638–641.  https://doi.org/10.1037/0033-2909.86.3.638
  * Simmons . (2011) simmons2011falseSimmons, JP., Nelson, LD. Simonsohn, U.  2011\.  False-positive Psychology: Undisclosed Flexibility in Data Collection and Analysis Allows Presenting Anything as Significant False-positive psychology: Undisclosed flexibility in data collection and analysis allows presenting anything as significant.  Psychological Science22111359–1366.  https://doi.org/10.1177/0956797611417632
  * Simonsohn . (2014) Simonsohn2014-cnSimonsohn, U., Nelson, LD. Simmons, JP.  2014\.  p-Curve and Effect Size: Correcting for Publication Bias Using Only Significant Results p-curve and effect size: Correcting for publication bias using only significant results.  Perspectives on Psychological Science96666–681.  https://doi.org/10.1177/1745691614553988
  * Stanley (2005) Stanley2005-ngStanley, TD.  2005\.  Beyond Publication Bias Beyond publication bias.  Journal of Economic Surveys193309–345.  https://doi.org/10.1111/j.0950-0804.2005.00250.x
  * Stanley Doucouliagos (2014) stanley2014metaStanley, TD. Doucouliagos, H.  2014\.  Meta-Regression Approximations to Reduce Publication Selection Bias Meta-regression approximations to reduce publication selection bias.  Research Synthesis Methods5160–78.  https://doi.org/10.1002/jrsm.1095
  * Sterling (1959) sterling1959publicationSterling, TD.  1959\.  Publication Decisions and their Possible Effects on Inferences Drawn from Tests of Significance—or Vice Versa Publication decisions and their possible effects on inferences drawn from tests of significance—or vice versa.  Journal of the American Statistical Association5428530–34.  https://doi.org/10.1080/01621459.1959.10501497
  * van Aert . (2016) Van_Aert2016-cuvan Aert, RCM., Wicherts, JM. van Assen, MALM.  2016\.  Conducting Meta-Analyses Based on p Values: Reservations and Recommendations for Applying p-Uniform and p-Curve Conducting meta-analyses based on p values: Reservations and recommendations for applying p-uniform and p-curve.  Perspectives on Psychological Science115713–729.  https://doi.org/10.1177/1745691616650874
  * van Assen . (2015) Van_Assen2015-qsvan Assen, MALM., van Aert, RCM. Wicherts, JM.  2015\.  Meta-analysis Using Effect Size Distributions of Only Statistically Significant Studies Meta-analysis using effect size distributions of only statistically significant studies.  Psychological Methods203293–309.  https://doi.org/10.1037/met0000025
  * Vevea Hedges (1995) Vevea1995-onVevea, JL. Hedges, LV.  1995\.  A General Linear Model for Estimating Effect Size in the Presence of Bublication Bias A general linear model for estimating effect size in the presence of bublication bias.  Psychometrika603419–435.  https://doi.org/10.1007/bf02294384
  * Vevea Woods (2005) Vevea2005-xpVevea, JL. Woods, CM.  2005\.  Publication Bias in Research Synthesis: Sensitivity Analysis Using A Priori Weight Functions Publication bias in research synthesis: Sensitivity analysis using a priori weight functions.  Psychological Methods104428–443.  https://doi.org/10.1037/1082-989X.10.4.428
  * von Neumann (1951) von1951variousvon Neumann, J.  1951\.  Various Techniques Used in Connection with Random Digits Various techniques used in connection with random digits.  Applied Math Series1236-385.


in Ref. [28]).

Let ${\mathcal{N}_{p}^{(\mathrm{dep})}}$ be the single-qubit depolarizing
channel of rate $p$, i.e.,
${\mathcal{N}_{p}^{(\mathrm{dep})}}(X)=p\frac{I}{2}\Tr(X)+(1-p)X$. Let
$\Phi:=\bigcirc_{i=1}^{L}(\mathcal{N}_{p}^{(\mathrm{dep})\otimes
n}\circ\mathcal{U}_{i})$ be a circuit of $L$ unitary layers interspersed by
local depolarizing noise. Then, for every state $\rho$, we have

$\displaystyle
D_{2}\left(\Phi(\rho)\bigg{\|}\frac{I}{2^{n}}\right)\leq(1-p)^{2L}n.$ (281)

Then the desired result follows by a direct application of the Pinsker’s
inequality.

######  Proposition 67 (Deviation from maximal mixedness).

Under the same assumptions of Lemma 66, we obtain

$\displaystyle\bigg{\|}\Phi(\rho)-\frac{I}{2^{n}}\bigg{\|}_{1}\leq\sqrt{2n}(1-p)^{L}.$
(282)

And therefore, for all states $\rho$ and $\sigma$, we have

$\displaystyle\norm*{\Phi(\rho)-\Phi(\sigma)}_{1}\leq 2\sqrt{2n}(1-p)^{L}.$
(283)

###### Proof.

Combining the Pinsker’s inequality with the monotonicity of the the quantum
Rényi divergence, we obtain

$\displaystyle\bigg{\|}\Phi(\rho)-\frac{I}{2^{n}}\bigg{\|}_{1}^{2}\leq
2D\bigg{(}\Phi(\rho)\bigg{\|}\frac{I}{2^{n}}\bigg{)}\leq
2D_{2}\bigg{(}\Phi(\rho)\bigg{\|}\frac{I}{2^{n}}\bigg{)}$ (284)

Hence Proposition 67 implies

$\displaystyle\bigg{\|}\Phi(\rho)-\frac{I}{2^{n}}\bigg{\|}_{1}\leq\sqrt{2n}(1-p)^{L}.$
(285)

Thus, a direct application of the triangle inequality yields the desired
result:

$\displaystyle\norm*{\Phi(\rho)-\Phi(\sigma)}_{1}\leq\bigg{\|}\Phi(\rho)-\frac{I}{2^{n}}\bigg{\|}_{1}+\bigg{\|}\Phi(\sigma)-\frac{I}{2^{n}}\bigg{\|}_{1}\leq
2\sqrt{2n}(1-p)^{L}.$ (286)

∎

We emphasize that the above result does not require randomness, unlike our
Theorem 21. On the other hand, Theorem 21 yields non trivial bound even at
constant depth, while the above depolarizing-noise result is vacuous at sub-
logarithmic depth. We remark that a non-vacuous bound for shallow depths could
be derived by means of the quantum Bretagnolle-Huber inequality (see, for
instance, Ref. [114], (Lemma B.1 in Ref. [115]) and references therein).

Informally, Proposition 67 says that the output of a noisy circuit becomes
computationally trivial at super-logarithmic depths, provided that the noise
is unital and primitive. This poses severe constraints on the capabilities of
noisy devices, as exemplified by the following result.

######  Corollary 68 (Exponential concentration, unital case).

Under the same assumptions of Lemma 66, for every state $\rho$ and for every
observable $O$, we obtain

$\displaystyle\left|\Tr[O\Phi(\rho)]-\frac{\Tr[O]}{2^{n}}\right|\leq\sqrt{2n}(1-p)^{L}\norm*{O}_{\infty}.$
(287)

###### Proof.

First, we notice that the LHS can be rearranged as follows:

$\displaystyle\left|\Tr[O\Phi(\rho)]-\frac{\Tr[O]}{2^{n}}\right|=\left|\Tr\left[O\left(\Phi(\rho)-\frac{1}{2^{n}}\right)\right]\right|.$
(288)

Hence, we obtain

$\displaystyle\left|\Tr\left[O\left(\Phi(\rho)-\frac{1}{2^{n}}\right)\right]\right|\leq\left\|\Phi(\rho)-\frac{I}{2^{n}}\right\|_{1}\norm*{O}_{\infty}\leq\sqrt{2n}(1-p)^{L}\norm*{O}_{\infty},$
(289)

where the first inequality follows from the Holder’s inequality, and the
second one is a consequence of Proposition 67. ∎

###  E.2 Numerical simulations

In this section we present numerical results that corroborate our analytical
results and explore regimes that go beyond the assumptions and results of our
theorems, such as the one of assuming the $2$-qubit gates to be distributed
according to a local $2$-design. We start by considering a brickwork
architecture as depicted in Fig. 2, where each 2-qubit gate takes the form

$\displaystyle
U_{i,i+1}\left(\theta_{1},\theta_{2},\theta_{3},\theta_{4}\right)\coloneqq\left(R_{Y}(\theta_{4})\otimes
R_{Y}(\theta_{3})\right)\mathrm{CNOT}_{i,i+1}\left(R_{X}(\theta_{2})\otimes
R_{X}(\theta_{1})\right).$ (290)

We consider the noise model given by the composition of amplitude damping and
depolarizing channels:

$\displaystyle\mathcal{N}_{(p,q)}^{(\mathrm{dep,amp})}\coloneqq\mathcal{N}_{p}^{(\mathrm{dep})}\circ\mathcal{N}_{q}^{(\mathrm{amp})},$
(291)

where two noise channels are defined respectively in Eq. (75) and Eq. (74)
with $p,q\in[0,1]$.

Moreover, we assume, in contrast to our circuit model in Fig. 2, that the
circuit ends with a layer of noise (instead of a layer of single-qubits
gates). Furthermore, we consider an expectation value with respect to the
observable $Z_{1}$.

Figure 7: Left: Variance of the partial derivative with respect to a gate in
the $L$-th layer of $5$ qubits one-dimensional quantum circuit of depth $10$.
We confirm that the expectation value is exponentially less ‘sensitive’ to the
gates in the first part of the circuit compared to the last layers. We set the
depolarizing and amplitude damping parameters to $p=q=0.2$. If we switch off
the amplitude damping ($q=0$), we observe that all the partial derivatives
exhibit similar scaling uniformly across all the layers. All the free
parameters are uniformly sampled in $[0,2\pi]$. Right: Variance of the partial
derivative with respect to a gate in the last layer of the circuit versus the
number of qubits. The depth of the circuit is chosen to be $2n$, where $n$ is
the number of qubits. In the case of both amplitude damping and depolarizing
noise with $p=q=0.2$, we observe a consistent behavior where the partial
derivative with respect to the last layer remains constant regardless of the
number of qubits. However, when $q=0$, an exponential decay is observed. In
both plots, the variance was estimated based on $100$ randomly sampled points
in the parameter space, and all the free parameters are uniformly sampled in
$[0,2\pi]$.

In Fig. 7 (Left), we can clearly observe that the partial derivatives taken at
the end of the circuit are significantly larger compared to those taken at the
beginning of the circuit. This confirms the exponential decay that we proved
in Theorem 50. It is noteworthy that if we were to deactivate the amplitude
damping component, we would observe an average exponential concentration in
all partial derivatives, regardless of the layer at which the derivative is
taken. This aligns with the findings of the study on (depolarizing-)noise-
induced barren plateaus [20] and with our Proposition 57. Furthermore, in our
experiments, we observe that the partial derivatives in the final layers
remain constant as the number of qubits increases, as depicted in Fig. 7
(Right). This observation corroborates the conclusions drawn in our Corollary
49 In particular, this implies that the $2$-norm of the gradient remains
constant on average with respect to the number of qubits. Another
consideration is whether our theorems heavily depend on the assumption that
the $2$-qubit gate is sampled by a $2$-design, and whether this assumption can
be relaxed. Our evidence suggests that we might relax such assumption: in our
numerical simulations, we observe a similar trend of what we proved even with
more structured ansatz like QAOA, the _quantum approximate optimization
algorithm_ [94]), as demonstrated in Fig. 8. This provides further support for
the notion that the assumption regarding the $2$-qubit gates being sampled by
a $2$-design is not crucial.

Figure 8: The same quantities as in the previous figure. However, in this
case, we employed a QAOA circuit ansatz. Notably, we observed that the results
do not undergo a qualitative change. This provides evidence that our findings
are not limited to unstructured ansatz composed of $2$-design $2$-qubits
gates. The choice of parameters remains consistent with those used in the
previous plots, with the only difference being the use of $X_{1}$ as an
observable instead of $Z_{1}$ (as the expectation value of $Z_{1}$ should
always be zero in the noiseless case for any fixed QAOA circuit, due to
symmetry reasons). The QAOA ansatz chosen has the form
$\prod^{D}_{i=1}\exp(-i\beta_{i}H_{x})\exp(-i\gamma_{i}H_{z})\ket{+}^{\otimes
n}$, where $H_{z}\coloneqq\sum^{n-1}_{i=1}Z_{i}Z_{i+1}$ and
$H_{x}\coloneqq\sum^{n}_{i=1}X_{i}$, where each $\gamma_{i}$ and $\beta_{i}$
are sampled uniformly randomly in $[0,2\pi]$.