Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
149,119 characters · 31 sections · 188 citation commands
Literature Review and Evidence Aggregation: a Toolkit for Applied Micro
\singlespacing \pagenumbering{arabic}
Consider an analyst interested in estimating or predicting the size of an effect -- the estimand $\theta_0$. She has a dataset with estimates from several prior studies of similar causal effects with standard errors. We provide a toolkit for three goals that she might have:
To develop this toolkit, we alternate between theoretical discussions of methods and practical case-studies. We rely on an empirical Bayes approach to provide a unified theoretical framework. This framework assumes that the estimands $\theta_i$ are drawn from some distribution $\mu$, while the observed estimates are conditionally normal, $\hat \theta_i \sim N(\theta_i, \sigma_i^2)$. If there is selectivity, then observability of study $i$ might depend on the $Z$-statistic $Z_i = \hat \theta_i / \sigma_i$. The distribution $\mu$ might be of independent interest, but more importantly it serves as a device to construct predictions for $\theta_0$ for some new context.
To illustrate these methods, we present both small case studies with three to five prior estimates (drawing on sager2025, blundell2025, bailey2025), and larger case studies with dozens of estimates. The larger case studies are drawn from several fields of applied economics, studying the effect of active labor market policy card2018works, unemployment insurance cohen2025disemployment, nudges dellavigna2022rctstoscale, and unconditional cash transfers crosta2024unconditional.
\paragraph{Question 1: Comparing estimates to the prior literature} This question appears in every applied micro paper that estimates a quantity also studied in prior work. The analyst begins by asking whether the estimates in the prior literature differ primarily because of sampling variation in $\hat\theta_i$ or because of between-study variation in underlying effects $\theta_i$. Between-study variation can be measured by comparing the dispersion of estimates $\hat{\theta}_i$ to the standard errors $\sigma_i$. Then, she can compute the precision-weighted mean of the prior estimates.
An estimate of the mean and distribution of the underlying effects is useful in several scenarios. Suppose an analyst has a new estimate which addresses a bias in the existing literature; she can test whether her new estimate could have emerged from the distribution of effects in the existing literature. Alternatively, suppose that the analyst has a new better-powered identification strategy and she wants to quantify the improvement in precision. Finally, suppose the prior literature contains meaningful between-study variation and the analyst has a new estimate to add to this literature. She can point to meaningful between-study variation as a motivation to produce estimates for additional contexts and to search for mechanisms that predict heterogeneity.
We illustrate that each of these three scenarios actually arises using papers published in American Economic Journal: Economic Policy (AEJ Policy) in 2025. In fact, these scenarios appear to emerge frequently. Of the 61 papers published in 2025 in AEJ Policy, 27 papers discuss prior estimates of similar targets from the literature. A crude lower bound is therefore that meta-analysis can be used---without the author finding any additional related studies---in nearly half of recently-published research in one prominent journal in applied micro.
In settings with dozens of estimates, an analyst can go further. If the underlying distribution of effects is normal, then the empirical Bayes leads to linear shrinkage estimators of effects. However, she can also use richer models which do not assume that the underlying effects are distributed normally. Using maximum likelihood estimation, they can estimate a $t$ distribution to account for fat tails, or even use a nonparametric distribution to account for skewness or multi-modality. Empirically, we find that we can reject the normal model in favor of a richer model in three out of our four larger case studies.
\paragraph{Question 2: Predicting effects for new contexts} The second question---what is the predicted effect $\theta_0$ in a new context---can be addressed by meta-regression or a generalization using Gaussian Process priors. The classic linear meta-regression approach to this question involves simply regressing estimates on covariates and using the coefficients from the regression for a new set of covariate values. Meta-regression's prediction for the new context can be rewritten as a weighted average of the prior estimates. Reweighting is attractive because it gives the analyst a transparent understanding of how the existing studies contribute to a prediction for the new $\theta_0$.
Gaussian Process priors provide a way to predict $\theta_0$ which maintains transparency through weighted averages, but does not require the restrictive linear additivity assumption of meta-regression. Instead, the researcher chooses a hyperparameter which corresponds to how much they are willing to import insights from existing studies in distinct contexts. When a large number of prior estimates are available for the same context, meta-regressions and Gaussian Process priors make similar predictions. However, when a researcher wants to extrapolate to a new context with no prior studies or only a few prior studies, the distinction between the methods becomes meaningful.\\
We illustrate the predictions from these methods using a dataset on the employment effects of active labor market policies from card2018works. Consider an analyst who seeks to predict the long-term effect of training for the long-term unemployed (LTU). This prediction exercise necessarily involves extrapolation, because the review does not contain any studies of this exact horizon by program by population cell. Fortunately, there are related studies of: (i) the long-term effects of training for a different population (unemployment benefit recipients) and (ii) the short-term effects of training for LTUs.
The analyst may rely on these studies interpreted through the lens of a meta-regression or a Gaussian Process prior. Meta-regression or a Gaussian Process prior with a long length scale---meaning that the analyst is willing to draw heavily on studies from contexts that are further away---deliver a precise posterior for the program's forecasted effect. The standard deviation for the posterior mean is 1 p.p. or less, which is quite small relative to the a posterior mean of roughly 10 p.p. In contrast, if the analyst uses a Gaussian Process prior with a short length scale, the distribution for the prediction is far more uncertain, with a posterior standard deviation of 7 p.p. We discuss how the researcher's preferred economic model of employment determination can be used to inform the statistical choice of the length scale hyperparameter.
\paragraph{Question 3: Detecting and correcting selectivity} Turning to the third question, we finally examine selectivity, where reporting and publication decisions can depend on an estimate's $Z$-statistic. We review two types of approaches to selectivity. We first discuss methods for testing for the presence of selectivity based on the distribution of $p$-values across studies. We then turn to methods which leverage the richer information contained in the joint distribution of estimates $\hat \theta_i$ and standard errors $\sigma_i$. Leveraging this richer information allows us to not only test for selectivity, but also to estimate its extent, and to obtain effect distributions that correct for selectivity. We provide graphical intuition for identification in this setting.
We find that selectivity is ubiquitous and that accounting for selectivity can lead to dramatically different estimates of parameters of interest in the context of our four larger case studies. First, insignificant findings are much less likely to be published: the probability of publishing an insignificant estimate is between 5% and 28% as compared to a positive significant estimate. Second, failing to correct for selectivity leads researchers to overstate the mean effect. The mean effect---after correcting for selectivity---is between 12% and 21% of the simple mean across the four applications we study. Ignoring selectivity leads to severely biased conclusions about the underlying effects.
\paragraph{Roadmap} The remainder of the article is structured as follows: Section (ref) introduces the formal setup, Section (ref) discusses aggregation, Section (ref) discusses aggregation with covariates, Section (ref) discusses selectivity, Section (ref) discusses how to combine selectivity and covariates, and Section 7 provides a cookbook for practitioners. Before introducing our formal setup, the remainder of this section provides a brief review of some of the relevant methodological literature.
\paragraph{Key theoretical references for our review} In the following, we draw on several literatures in statistical theory that provide a foundation for the more specific questions encountered in evidence aggregation. The empirical Bayes (or random-effects) framework that organizes our discussion throughout originates with robbins1956empirical, stein1981estimation, and Morris1983; We rely on both parametric and non-parametric implementations: kernel deconvolution meister2009, convex programming for non-parametric maximum likelihood of the prior koenker2014convex, Gaussian process priors for flexible extrapolation across covariates williams2006gaussian, and Tweedie-type shrinkage of individual estimates efron2011tweedie. An approach closely related to empirical Bayes imposes priors on hyper-parameters, yielding hierarchical Bayes estimators, cf. gelman2014bayesian; we discuss this approach in (ref).
A central focus of our review is selective reporting and publication of empirical findings. Concerns about selectivity in published economics research date back at least to stanley1989metaregression. Bias tests via meta-regression were introduced by egger1997bias and refined by stanley2017finding. brodeur2016star document bunching of reported test statistics just above conventional significance thresholds, an approach extended to a broader cross-section of identification strategies by brodeur2020methods; ioannidis2017power highlight the prevalence of severely underpowered findings and propose to focus on highly-powered studies. elliott2022detecting provide non-parametric tests for selectivity based on violations of monotonicity of the p-curve. Likelihood-based methods that estimate selection from the joint distribution of estimates and standard errors are developed in publicationbias2019; we build on these in Section (ref). christensen2018transparency survey the broader research-transparency movement concerned with selective reporting.
\paragraph{Meta-analysis in economics} Quantitative research synthesis has a long tradition in economics, with early applications such as CardKrueger1995's synthesis of the minimum-wage literature. Methodological tools for meta-regression and the detection of publication selection are developed by Stanley2008 and stanley2014meta, and surveyed in the practitioner's guide of irsova2023meta. A more recent wave exploits standardized designs and larger samples of studies, including the comparison of academic and at-scale nudge trials in dellavigna2022rctstoscale, the replication exercise of camerer2016evaluating, and the generalizability analysis of vivalt2019much.
\paragraph{Approaches beyond the scope of the present review} There are a number of exciting and insightful developments in the recent literature on evidence aggregation which are beyond the scope of the present review:
An increasing set of studies, in particular in development and experimental economics, make micro-data available. In settings where this is the case, meta-studies can directly analyze the pooled micro-data, yielding insights beyond those available from aggregate statistics; see, e.g., bandiera2021social, meager2022aggregating, kremer2023cleanwater, lund2024mentalhealth, mullins2025welfare.
Meta-regressions, as discussed in (ref), are often used to explore possible drivers of effect heterogeneity. Parts of the literature use priors for regression coefficients putting a point-mass at 0. This leads to Bayesian model averaging approaches irsova2023meta. In our review (and in line with the recommendations of gelman2014bayesian), we instead focus on continuous priors that do not assume true sparsity in a correctly specified linear model.
A central focus of our review is prediction of effects $\theta_i$ for specific instances. A decision-theoretic foundation for the proposed approaches can be found in the theory of compound decision problems zhang2003compound. An active literature in theoretical econometrics explores the broader decision-theoretic foundations and ramifications of evidence aggregation, e.g.,\ manski2020toward, christensen2026optimal, ishihara2024evidence, and montielolea2023decision.
Throughout, our review is based on the following framework. This framework allows us to provide a unified discussion of the issues involved in meta-studies, including information aggregation and extrapolation, as well as p-hacking and publication bias.
Suppose we are interested in a set of empirical studies indexed by $i \in 1,\ldots,n$. Corresponding to each of these studies there is an (unknown) estimand $\theta_i$. This might for example be the true average treatment effect, if $i$ is an experimental study, or it might be the elasticity of some behavioral relationship. Each study reports an estimate $\hat \theta_i$ of $\theta_{i}$. Each study is furthermore characterized by study-level covariates $X_i$, which might include features such as sampling frame and site of the study, study design and identifying variation, treatment or policy characteristics, etc. We will use $\theta$ to denote $(\theta)_{i\in1, \ldots, n}$ and use $\hat\theta, X$ analogously.
Not all studies are available to us when we conduct a meta-analysis. Some studies are never published, either because of researcher decisions (possibly reflecting issues such as p-hacking or specification searching), or because of reviewer and editor decisions (which might lead to publication bias). Some studies might also simply be excluded from the sampling frame of a meta-study, for instance based on the journal or year they were published in. We indicate observability and inclusion of a study in a meta-analysis by $D_i \in \{0,1\}$, where $D_{i}=1$ for studies that are included.
We assume that the sampling distribution of $\hat{\theta}_{i}$ is given by $$\hat \theta_i |\theta_i, \sigma_i^2,X_{i} \sim N(\theta_i, \sigma_i^2),$$ where the standard error $\sigma_{i}$ is known to us whenever a study is published. This assumption implies that $\hat{\theta}_{i}$ is unbiased for $\theta_{i}$; alternatively, we can think of this assumption as simply defining the estimand $\theta_{i}$ as the expectation of $\hat{\theta}_{i}$. This assumption also implies that $\hat{\theta}_{i}$ is normally distributed. Normality is typically justified by asymptotic arguments based on the central limit theorem. The same asymptotic arguments also imply that we can treat $\sigma_i^2$ as known.
Under normality, it is useful to define the $Z$-statistic $Z_{i}$, the corresponding normalized parameter $\omega_{i}$, and the p-value $P_{i}$ for a one-sided test of the null hypothesis that $\theta_{i}=0$, via $$
$$ where $\Phi$ is the cumulative distribution function of the standard normal distribution.
We furthermore assume that the estimands $\theta_{i}$ are themselves sampled from some unknown distribution $\mu$ across studies, $$ \theta_{i}\sim^{iid} \mu. $$ The distribution $\mu$ is the population counterpart to the sample distribution $\{\theta_{1}, \dots, \theta_{n}\}$. The sample distribution is well-defined, regardless of the substantive relationship between the estimands of different studies $i$. We discuss the interpretation of the distribution $\mu$ in more detail below in Section (ref).
We will at times rely on the stronger assumption $\theta_{i} |\sigma_i \sim^{iid} \mu$, so that the estimands vary independently of estimator precision across studies. For the most part, this assumption only serves as a convenient simplification for the construction of estimators. The one point where this assumption is truly required is identification of selection models using meta-studies, discussed in (ref).
The observability $D_{i}$ of studies, which reflects reporting and publication decisions, can in principle depend on $\hat{\theta}_{i}, \sigma_{i}$, and $X_{i}$ in complicated ways. We will however assume that the probability of being observable only depends on the $Z$-statistic $Z_{i}$, $$P(D_i =1 | \hat \theta_i, \sigma_i) = \bar{d}(Z_{i}).$$ This is a restrictive assumption, which is maintained in much of the literature on p-hacking and publication bias (cf. Section (ref)).
We start our methodological review by considering methods for evidence aggregation. For clarity of exposition, throughout this section we assume that there is no selective reporting, so that $\bar{d} \equiv 1$. The natural framework for discussing evidence aggregation is the empirical Bayes framework, which was first introduced by robbins1956empirical. The approach, as applied to our setting, can be summarized as follows:
Many meta-studies stop at step 3 of this recipe, characterizing properties such as the expectation $E[\theta_{i}]$ (the average effect), the variance $Var(\theta_{i})$ (effect heterogeneity), or the covariance $Cov(\theta_{i}, X_{i})$ (predictors of effect heterogeneity). Step 4 in this recipe is, however, necessary whenever the goal is to solve a specific policy problem, where policies are chosen for a specific site $i$, rather than across all possible sites.
Viewed through the lens of these four steps, the distribution $\mu$ is best thought of as a conceptual device. It is an input to the construction of empirical Bayes estimators, decision procedures in compound decision problems, and as a device to discuss questions of external validity. However, the justification of these procedures does not necessarily rely on the family of priors. In addition, we note that empirical Bayes has provably good properties for such compound decision problems, conditional on the $\theta_{i}$; see for instance stein1981estimation.\footnote{Compound decision problems involve repeated decisions across the instances $i$, and loss that is averaged across these instances. A leading example is estimation of the vector $(\theta_1, \dots, \theta_n)$, with mean squared error $\tfrac1n \sum_i (\hat \theta_i-\theta_i)^2$.} For a general review, see zhang2003compound. The primary object of interest for both analysts and policymakers will typically be the effect $\theta_0$ for new instances, rather than the distribution $\mu$ itself.
We start by considering the parametric empirical Bayes approach absent covariates, with a normal family of priors (settings with covariates are discussed in (ref) below). In doing so, we follow Morris1983; see also rubin1981estimation, and dersimonian1986meta. We next generalize by allowing that $\mu$ has heavy tails, considering a family of $t$-distributions. We can finally leave the distribution of $\theta_{i}$ fully unrestricted.
Let us first consider the simplest empirical Bayes setting, following Morris1983, where we ignore covariates and assume that the effects $\theta_i$ are normally distributed across studies, independently of $\sigma_i$, so that $\mu$ is given by $$\theta_i | \sigma_{i} \sim N(\bar{\theta}, \tau^2).$$ Under this assumption, the marginal distribution of $\hat{\theta}_{i}$, integrating out the latent $\theta_{i}$, is given by $$\hat{\theta}_{i}|\sigma_{i} \sim N(\bar{\theta}, \tau^{2} + \sigma_{i}^{2}).$$ A simple method-of-moments estimator of $\bar{\theta}$ is given by $\hat{\bar{\theta}} = \tfrac{1}{n}\sum_{i} \hat{\theta}_{i}$. A simple estimator of $\tau^{2}$ is given by
The efficiency of the estimator of $\bar{\theta}$ could then be further improved by using a precision-weighted average,
Alternatively, both $\bar{\theta}$ and $\tau^{2}$ can be estimated simultaneously using maximum likelihood. In the absence of any further information about a new instance $\theta_0$, we can think of $\hat{\bar{\theta}}$ as a predictor for $\theta_0$.
To estimate $\theta_{i}$ for individual instances $i$ where we do have an estimator $\hat{\theta}_{i}$, we can then apply step 4 of the empirical Bayes recipe, where
Substituting estimates for $\bar{\theta}$ and $\tau^{2}$ into this expression yields the empirical Bayes estimator of $\theta_i$. It is interesting to note that in the special case where $\sigma^{2}_{i}$ is constant across $i$, the empirical Bayes approach essentially recovers the James-Stein shrinkage estimator (up to a small degrees of freedom correction). This shrinkage estimator is known to uniformly dominate $(\hat{\theta}_i)_{i=1}^n$ as estimator of the vector of $\theta_i$'s; see stein1981estimation for a proof.
The normality assumption $\theta_i \sim N(\bar{\theta}, \tau^2)$ underlying equation (ref) is restrictive: in many applications the latent distribution $\mu$ has heavy tails, is skewed, or is multi-modal. We can generalize our baseline model to relax this assumption.
First, one can take $\mu$ to be a (possibly scaled and shifted) $t$-distribution with $\nu$ degrees of freedom,
which allows for heavier tails than the normal; smaller $\nu$ implies heavier tails.
Second, one can leave $\mu$ entirely unrestricted. In this case, $\mu$ can be estimated by nonparametric maximum likelihood koenker2014convex or by deconvolution meister2009, and the posterior mean of $\theta_i$ admits a model-free expression known as Tweedie's formula efron2011tweedie:
where $f(\hat{\theta}_i\,|\,\sigma_i)$ is the marginal density of $\hat{\theta}_i$ given $\sigma_i$. Intuitively, the posterior mean adjusts the raw estimate $\hat{\theta}_i$ in the direction of higher marginal density: the magnitude and direction of shrinkage are determined by the local shape of $f$ rather than by a global normality assumption. We illustrate the empirical relevance of this generalization, which allows us to appropriately adjust for heavy tails, skewness, and multi-modality, in Section (ref). The derivation of Tweedie's formula and its connection to deconvolution are given in Appendix (ref).
The primary focus of many applied micro papers is to credibly estimate a quantity of interest and then to interpret the new estimate in part by comparing it to prior work. To illustrate how meta-analysis might be able to add additional insights to this type of exercise, we start from the studies which a paper already cites as estimating a similar object and ask what additional insights can be extracted based on applying the equations from Section (ref).
We distinguish between three uses of meta-analysis. In the first example, an analyst has what they believe to be a superior methodology and wants to know if their estimate is meaningfully different from the prior literature. In the second example, an analyst has a new estimate and they want to construct the most likely estimate for the true effect, taking into account the distribution of effects in the prior literature. In the third example, an analyst has a new estimate and they want to quantify the improvement in precision.
\paragraph{Case Study \#1: Correcting Bias} sager2025 explore different strategies to estimate the effect of pollution on housing prices, with an eye towards correcting bias in the methods used in the prior literature. They describe their paper as contributing three insights to the literature on pollution damages. After describing the first two insights, they write
A meta-analytic framework offers additional lessons which complement the discussion in sager2025. With the estimates from Table (ref), one could therefore have expanded on the “third insight” discussion in that paper, adding
\paragraph{Case Study \#2: Best Estimate of True Effect} blundell2025 is interested in estimating the effect of pay transparency laws on the gender pay gap. They study a 2018 UK law that requires firms with over 250 employees to publicly disclose “gender equality indicators.” When discussing their paper in the context of the broader literature, they mention three other papers that attempt to estimate a similar quantity, writing:
Meta-analysis yields additional insights which help to put this paper's findings in context. Unlike in sager2025, blundell2025 believes the prior literature to be relevant to their context, so this new estimate is the most efficient estimate available. With the estimates from Table (ref) one could therefore have expanded on the discussion in that paper, adding
The three examples in this section are drawn from articles published in American Economic Journal: Economic Policy in 2025. To assess the generalizability of the examples listed here, we reviewed all 61 articles published in the journal in 2025. In essentially all of these articles, the primary object of interest is a causal effect; however, many of the articles do not directly compare their estimates to prior papers or only compare their estimates to a single prior paper. Table (ref) lists 27 papers which compare their new estimates to at least two prior estimates from the literature. A crude lower bound is therefore that meta-analysis can be used---without identifying any additional related studies---in almost half of high-quality research in one journal in applied micro.
\paragraph{Case Study \#3: Improving Precision}
bailey2025 is interested in the effect of paid family leave on mothers' earnings and employment. They use a regression discontinuity design in California to estimate the impact of this policy. Describing the prior literature, they write
Table (ref) conducts a meta-analysis for the effect of offering paid family leave on earnings. With such a meta-analysis in hand, the discussion might have continued
The simple meta-analyses done in this section are likely to generalize to other research in applied micro. All three examples are drawn from articles published in American Economic Journal: Economic Policy in 2025. To assess the generalizability of the examples listed here, we reviewed all 61 articles published in the journal in 2025. In essentially all of these articles, the primary object of interest is a causal effect; however, many of the articles do not directly compare their estimates to prior papers or only compare their estimates to a single prior paper.
Appendix (ref) lists 27 papers which compare their new estimates to at least two prior estimates from the literature. A crude lower bound is therefore that meta-analysis can be used---without identifying any additional related studies---in almost half of high-quality research in one journal in applied micro.
\paragraph{Applications} In the remainder of the paper, we re-analyze datasets from four existing meta-analyses: card2018works, cohen2025disemployment, crosta2024unconditional, and dellavigna2022rctstoscale. These meta-analyses cover four different subfields of economics -- public, development, behavioral, and labor -- which allows us to emphasize considerations specific to different context. Unlike the “small” meta-analyses in the previous section, each of these meta-analyses is based on a relatively large number of studies and uses covariates to predict treatment effects. The larger $n$ and the covariates allow us to use richer methodologies to extract additional insights.
cohen2025disemployment analyzes the effect of unemployment benefit generosity on unemployment duration using 93 independent estimates from 57 studies, most of which are based on natural experiments. There are no randomized controlled trial (RCT) experiments in the sample, and a few older studies rely on a selection-on-observables assumption rather than a natural experiment. The paper finds that insignificant results are about 10% as likely to be published as significant results. The paper also finds that the baseline replacement rate for unemployment benefits is a meaningful predictor of the benefit duration elasticity, a point to which we return in Section (ref).
card2018works analyzes the effect of active labor market policies (ALMPs) on employment. About one-fifth of the studies the paper considers use an RCT design, and more than half of the studies occurred in Europe. The paper has an extensive heterogeneity analysis by covariates, which we revisit in Section (ref). Standard errors are missing from many of the studies in the original meta-analysis. We drop those studies for a final sample of 169 estimates from 45 studies.
dellavigna2022rctstoscale analyzes the effect of nudges on a wide range of outcomes. The paper considers the effect of RCTs from two different samples: 74 estimates that were published in academic journals and 241 estimates from studies conducted by two government-run “nudge units.” dellavigna2022rctstoscale find substantial publication bias in the sample of academic studies, relative to the nudge unit sample. We will combine both the nudge unit and academic journal estimates. In the selectivity section, we will assume no selectivity, $\bar d \equiv 1$, for estimates from the nudge unit subsample.
crosta2024unconditional analyzes the effect of unconditional cash transfer programs in low- and middle-income countries on several different outcomes. The paper exclusively considers RCT studies. crosta2024unconditional concludes, based on 75 estimates, that for each additional \$100 of transfer, recipients have monthly consumption that is \$3.30 higher.
\paragraph{Evidence aggregation} (ref) illustrates four lessons about evidence aggregation. First, for all four applications, there is substantial dispersion in estimated treatment effects beyond what would be expected based only on standard errors. This is reflected in the fact that we estimate $\tau^2 > 0$ both when using a normal distribution or a $t$ distribution for the latent distribution of estimates.
Second, using the MLE normal (or equation (ref)) provides a principled way to reduce weight on noisy outlier estimates. This is most easily seen in the context of crosta2024unconditional, where the simple mean and MLE normal deliver meaningfully different conclusions. The simple mean implies that a \$1 transfer raises consumption by \$0.13 while MLE normal implies that \$1 transfer raises consumption by \$0.03. This pattern arises because crosta2024unconditional is a setting where a few estimates have unusually large treatment effects and standard errors. The simple mean weights these outliers equally to the other estimates. Indeed, when the 10 percent of estimates with the largest standard errors are excluded, the method of moments estimator shrinks to \$0.04, in line with the MLE normal estimate.
Third, MLE estimates of the mean can be much more precise than the method of moments estimates. The precision advantages of the MLE are greatest when there is (i) variation in sampling uncertainty ($\sigma_i^2$) and (ii) between-study variation is much smaller than sampling uncertainty ($\tau^2<<\sigma_i^2$). This is illustrated in the estimates for crosta2024unconditional, which meets both of these conditions. At the opposite extreme, when $\tau^2 >> \sigma_i^2$, as in dellavigna2022rctstoscale, the weights on each study will be similar, leading MLE normal to yield a similar standard error to the method of moments. Standard errors for MLE estimators are calculated using a sandwich estimator. We include standard errors for $\hat\tau$ in the “large” meta-analyses in (ref) but not in the “small” meta-analyses in Tables (ref) and (ref) since these standard errors are imprecisely estimated in small samples.
Fourth, a non-normal model often provides a better fit to the data. For three of the applications (cohen2025disemployment, card2018works, and dellavigna2022rctstoscale), we reject that the latent distribution is normal. This is reflected in the finite degrees of freedom parameter $df$. In these three applications, assuming that the latent distribution is $t$ leads to a much lower estimate for the latent median $\bar{\theta}$ (as compared to assuming that the latent distribution is normal). This finding will serve as a justification for why we prefer to use a $t$ distribution in the context of analyzing selectivity in Section (ref). For crosta2024unconditional, assuming the latent distribution is $t$ leads to pathological results since the normal model is the best fit.
Until now, we have studied models which assume that the latent distribution of effects can be modeled using a normal or a $t$ distribution. When the latent distribution is in fact skewed or multi-modal, these parametric assumptions will lead to biased conclusions. A more flexible alternative---the nonparametric MLE koenker2014convex---combined with Tweedie's formula (equation (ref)) allow the shape of the posterior mean to adapt to the data.
(ref), Panel A, shows shrinkage curves under three simulated latent distributions that violate normality: a bimodal mixture of two normal distributions, a heavy-tailed $t_{3}$, and a right-skewed (centered) $\chi^{2}_{3}$. In each case the NPMLE curve departs sharply from the normal benchmark in the tails. In the bimodal case it is strongly non-monotone in slope, pulling observations toward the closer mode rather than a single global average. The $t$-model captures heavy tails but cannot recover skewness or multimodality. In the skewness case, right-tail estimates far away from zero are shrunk less than closer estimates. Panel B shows the analogous comparison for the meta-analysis of unemployment elasticities from cohen2025disemployment: the marginal density of $\hat\theta_i$ is visibly right-skewed with a heavier right tail than a matching normal, and the NPMLE shrinkage curve is correspondingly asymmetric---it shrinks negative estimates aggressively but flattens out for large positive $\hat\theta_i$ rather than pulling them back to the global mean.
The choice of latent distribution is most consequential for study-level estimates of individual $\theta_{i}$---ranking treatments, identifying “best” effects, or cutoff-based decisions. For aggregate quantities like $\bar{\theta}$, the differences across models are typically modest. As a practical guideline: when $n\gtrsim 50$ and the empirical distribution of $\hat \theta_i$ shows visible deviations from normality, the NPMLE (or at least a $t$-model) should be preferred; with small $n$, parametric models might remain more reliable.
Thus far, we have studied what a meta-analysis can discover based on estimates and standard errors alone. Dispersion in the distribution of effects across studies ($\tau^{2}>0$) naturally motivates a search for predictors of study-level heterogeneity. In each of our four applications, readers and researchers will naturally wonder why treatments are more effective in some instances than others. In this section, we consider what more can be learned when study-level covariates are available.
The basic empirical Bayes approach can be generalized by allowing for covariates, which enter as linear predictors of $\theta_{i}$; see for example chetty2014land. A further generalization, using tools that might be less familiar to readers in economics, allows for a non-linear relationship between the covariates $X_{i}$ and the conditional mean of $\theta_{i}$ given $X_i$, and does not restrict the functional form of this relationship. This can be implemented using Gaussian process priors williams2006gaussian, which allow for closed form predictions, despite the unrestricted functional form.
This approach is particularly useful when we want to think about the external validity of a meta-study, and the extrapolation of estimates to new problem instances. This approach allows for an assessment of (posterior) uncertainty, when extrapolating outside the support of previously encountered instances. By contrast, extrapolations based on (linear) functional forms might give a false sense of certainty.
To make predictions of $\theta_0$ for new instances, one needs to learn the relation of effect heterogeneity to contextual covariates $X_i$, to make predictions for new instances characterized by $X_0$. This motivates meta-regressions of the form
as estimated by stanley1989metaregression and by many authors since. We emphasize that the variation in $X_i$ is not exogenous and so this regression identifies a correlation between the covariate and the treatment effect rather than a causal relationship.
To rationalize such meta-regressions, and generalizing the basic parametric empirical Bayes model, assume that
Here the coefficient vector $\beta$ takes the role previously held by $\bar{\theta}$. For this model, we then get $$\hat{\theta}_{i} |X_{i} = x, \sigma_{i} \sim N(x \cdot \beta, \tau^{2} + \sigma_{i}^{2}).$$
Generalizing the method-of-moments approach of Section (ref), $\beta$ can be estimated using either OLS, or a weighted least squares (WLS) regression where the weights are $\frac{1}{\sigma_i^2}$. As a further generalization, we can run a regression of $\hat{\theta}_i$ on $x$, and $\tau^{2}$ can then be estimated as the average of $e_i^{2} - \hat{\sigma}_i^{2}$, where $e_i$ are the estimated residuals of the WLS regression. We could then re-run the regression, weighting each observation by the inverse of its variance, $\frac{1}{\tau^2 + \sigma_i^2}$, which yields a weighted least squares random effects (WLS RE) predictor.
Alternatively we could estimate the hyper-parameters $\tau^2$ and $\beta$ jointly by maximizing the marginal likelihood, as in the empirical Bayes approach. Or we could impose a prior on both $\tau^2$ and $\beta$ to obtain a hierarchical Bayes model, as discussed further below.
It is common to think of predictions based on regressions as a two-step process: First, the regression coefficients $\beta$ are estimated. Then the effect $\theta_0$ for a new instance is predicted as $X_0 \cdot \widehat{\beta}$. There is an alternative, numerically equivalent way to obtain the same result: we can leap-frog estimation of $\beta$ and directly predict $\theta_0$.
This alternative perspective, which we will explain next, has two advantages for our purposes. First, it makes explicit that predictions are weighted averages of the available estimates in a meta-study, and clarifies the weights used. Second, it lends itself naturally to generalizations from the linear model to non-linear (but smooth) relationships between $\theta$ and $X$. Using such more general models in turn is helpful because it avoids undue reliance on functional form for extrapolation.
\paragraph{Prediction in the linear model (Ridge)}
In this spirit, let us now characterize the posterior for $\theta_{0}$ as a function of $ (\hat{\theta}_{1}, \dots, \hat{\theta}_{n})$. In addition to the linear model of Equation (ref), assume that $\tau^2$ is known (possibly from first-stage maximum likelihood estimation), and suppose that we have a prior for the coefficient vector $\beta$ of the form $$ \beta \sim N(0, \rho^2 \cdot I), $$ where $\rho^2$ is known. Then the prior distribution of $(\theta_{0}, \hat{\theta}_{1}, \dots, \hat{\theta}_{n})$ is jointly normal conditional on $\tau^2$ and conditional on the observed covariates $(X_0, X_1, \dots, X_n)$; we will leave the conditioning implicit in our notation. The prior moments of this joint normal distribution are given as follows, for $i,j=0,1,\dots,n$: $$
$$ where $\rho^{2} + \tau^{2}$ governs the prior variance of $\theta_i$. The ratio $\rho^{2} / ( \rho^{2} + \tau^{2})$ determines the expected share of variation of $\theta_i$ that is explained by $X_i$ (analogous to an $R^2$ statistic). We collect these terms in the $n$-vector of prior covariances $c$ and the $n\times n$ prior variance matrix of the estimands $C$,
to get
Given joint normality of $(\theta_{0}, \hat{\theta}_{1}, \dots, \hat{\theta}_{n})$, the posterior mean of $\theta_0$ is equal to the posterior best linear predictor,
For intuition, note that this formula for the best linear predictor is the multivariate counterpart of the familiar formula for the univariate OLS slope, which is equal to the covariance of prediction target and predictor, divided by predictor variance. The posterior distribution of $\theta_{0}$ is thus given by
Note that the posterior mean of $\theta_{0}$ is a weighted average of the $\hat{\theta}_{i}$, with larger weights put on studies $i$ that have covariate values $X_{i}$ closer to $X_{0}$ (so that $C_{0,i}= \rho^2\cdot X_0 X_i'$ is larger), and on studies that have smaller values of $\sigma_{i}$. The covariances $C_{0,i}$ measure how similar (close) instance $0$ is to instance $i$. The coefficients $c \cdot (C + diag(\sigma_{i}^{2}))^{-1}$ then give the optimal weights for a posterior prediction. The variance $C_{0,0} - c \cdot (C + diag(\sigma_{i}^{2}))^{-1} \cdot c'$ measures the posterior uncertainty for $\theta_{0}$.
\paragraph{Nonlinear generalization (Gaussian Process)}
The linear model that we just discussed specifies the prior co-variance of $\theta_i$ and $\theta_j$ as $C_{ij}=K(X_i, X_j)$ for $i\neq j$ and $C_{ii}=K(X_i, X_i) + \tau^2$, where $K(x,x') = \rho^2\cdot x x'$ . The function $K(x,x')$ is called a covariance kernel. As we saw in equation (ref), the predicted value of $\theta_0$ can be calculated using just this covariance kernel, without any explicit reference to the regression slopes $\beta$.
We have assumed thus far that the prediction function $\bar{\theta}(x) = E[\theta_i | X_i=x]$ is linear; this is reflected in the prior covariances $K(x,x') = Cov(\bar{\theta}(x), \bar{\theta}(x'))$. This approach can be generalized beyond linear models. If we assume that the $\theta_i$ are jointly normal (conditional on covariates), with prior covariances specified by a general covariance kernel $K(X_i, X_j)$, then we obtain a so-called Gaussian process prior for $\bar{\theta}(x)$ williams2006gaussian. Denoting again $C_{ij}=K(X_i, X_j)$ for $i\neq j$ and $C_{ii}=K(X_i, X_i) + \tau^2$, the exact same expression as before (Equation (ref)) describes the joint prior distribution of $\theta_0$ and the observed estimates $(\hat{\theta}_{1}, \dots, \hat{\theta}_{n})$. Correspondingly, the same expression (Equation (ref)) continues to describe the posterior distribution of $\theta_0$ given the observed estimates.
The linear model specifies $K(x,x') = \rho^2\cdot x x'$. A popular alternative covariance kernel, which does not impose linearity, is the squared exponential kernel,
The “length scale” $l^{2}$ governs the smoothness of the relationship between covariates $X$ and effects $\theta$, where larger $l$ imply smoother functions. This length scale thus place a role analogous to the bandwidth in kernel regressions. These hyper-parameters $\rho^{2}, \tau^2$ and $l^2$ could be chosen a priori, for a fully Bayesian specification, or in a data-dependent manner, as in the empirical Bayes approach.
There are many other kernels which have been used in the literature. Spline regression, for example, corresponds to a Gaussian Process prior where $\bar{\theta}(\cdot)$ follows an integrated Brownian motion. While the squared exponential kernel leads to infinitely differentiable predictions, spline regressions are only twice differentiable, in general.
The empirical Bayes approach estimates hyper-parameters, such as $\tau^2$ and $\beta$ for the linear model, by looking at the data. Alternatively, one might put a prior on hyper-parameters (i.e. $\beta \sim \pi_\beta, \quad \tau \sim \pi_\tau$) and obtain a posterior distribution for hyper-parameters and study-level estimates in one step. The posterior based on the meta-regression model from (ref) can be written as
$$P(\theta, \beta, \tau |\hat\theta, \sigma_i, X) \propto \pi_\beta\pi_{\tau}\prod \frac{1}{\sigma_i}\varphi\left(\frac{\hat\theta_i - \theta_i}{\sigma_i}\right)\frac{1}{\tau}\varphi\left(\frac{\theta_i - X_i\beta}{\tau}\right).$$ This yields hierarchical Bayes estimators (see e.g. meager2019understanding).
If we put a normal prior $\beta \sim N(0, \rho^2 I)$ on the coefficients and a normal prior with very large variance on $\tau$, we recover the same estimated $\beta$ as the linear GP model (Ridge) discussed above. In particular $$ E[\beta | \hat{\theta}_1,\ldots, \hat{\theta}_n, \tau^2] = \operatorname*{argmin\;}_b \sum_i \frac{(\hat{\theta}_i-X_i\cdot b)^{2}}{ \tau^{2}+\sigma_i^{2} } + \frac{\|b\|^2}{\rho^2} . $$ We use Hierarchical Bayes in (ref) and (ref).\footnote{We also use a variant of Hierarchical Bayes in Section (ref) when we do estimation with Gaussian Process priors. In that case, we choose the hyper-parameters manually. This can be thought of as a prior which is a point mass (rather than a distribution).} One advantage of the hierarchical Bayes approach is that it facilitates uncertainty quantification for predictions, taking into account hyper-parameter uncertainty. This matters especially in the context of smaller meta-studies, where hyper-parameter uncertainty is non-negligible.\footnote{An alternative route to uncertainty quantification, which stays within the empirical Bayes framework rather than imposing priors on hyper-parameters, is provided by the robust empirical Bayes confidence intervals of armstrong2022robust.}
A leading downstream use of the apparatus developed above is the choice of binary implementation decisions $A_{i}\in \{0,1\}$ -- for example, whether to scale up a treatment evaluated in study $i$, or whether to implement a treatment in a new instance characterized by covariates $X_{0}$. If $\theta_{i}$ captures all relevant costs and benefits of the treatment, so that welfare is $U(A_{i},\theta_{i}) = A_{i}\cdot \theta_{i}$, a hypothetical decision-maker who observed $\theta_{i}$ might choose $A_{i} = \mathbf 1(\theta_{i}\geq 0)$. In practice we observe only $(\hat\theta_{i}, \sigma_{i}, X_{i})$, and the optimal Bayesian decision becomes
where the posterior mean is computed using the empirical Bayes or Gaussian process expressions above. Note that the sign of this posterior mean may differ from that of $\hat\theta_{i}$ when $\hat\theta_{i}$ is small relative to its standard error or inconsistent with the effects predicted by the covariates. For a new problem instance for which $\hat\theta_0$ is not observed (but covariates $X_{0}$ are observed), the analogous rule is $A_{0} = \mathbf 1\bigl( E[\theta_{0}\mid X_{0}, \hat{\theta}_1,\ldots, \hat{\theta}_n] \geq 0\bigr)$, using the Gaussian process posterior mean for $\theta_{0}$ from Equation (ref). We discuss other implementation decisions in more detail in Appendix (ref).
card2018works studies how the effectiveness of active labor market programs (ALMPs) varies with the program type (training, job search assistance, subsidized private sector employment, subsidized public sector employment), with the gender targeted (mixed-gender, female, male), with the time horizon (short-term, medium-term, long-term), and with the population targeted (UI recipients, long-term unemployed, or disadvantaged). The paper considers models where the first group in each parenthetical constitutes the reference group and asks how changes in program type, gender targeted, time horizon, and population targeted affect the probability of employment. The models used in the paper give each estimate equal weight when assessing the role of covariates.\footnote{This estimation approach is quite common. card2018works and dellavigna2022rctstoscale do this via Ordinary Least Squares, while cohen2025disemployment does this via Bayesian Model Averaging.} We capture several of the research questions examined in their paper by estimating
Figure (ref) shows how the effectiveness of ALMPs varies with these characteristics, comparing estimates from four types of models. The first regression weights every estimate equally and includes studies with missing standard errors, as in the original card2018works analysis.\footnote{In our analysis below, we treat the $\theta_i$ as independent draws. As noted in Section (ref), in the studies reviewed by card2018works, there are often multiple estimates $i$ per study $j$. Generalizing our approach, these might be modeled as $\theta_{ij} = \bar \theta_i + \eta_{ij}$, with each of the $\bar \theta_i$ and $\eta_{ij}$ independent. We do not pursue this here.} Using this model, the paper concludes that programs which target the long-term unemployed are more effective at raising employment. We replicate this finding (shown in purple diamonds): targeting the long-term unemployed raises employment effects by about 7 percentage points. The second limits the sample to studies with non-missing standard errors; the finding holds in this subsample as well.
Applying the frameworks from Section (ref), we explore alternative models that account for differences in precision across estimates and allow for a latent distribution of treatment effects. The third (“WLS”) weights each estimate by the inverse of its variance, $w_i = 1/\sigma_i^2$. The fourth (“WLS with random effects”) additionally allows for a latent distribution of treatment effects, weighting each estimate by $w_i = 1/(\tau^2 + \sigma_i^2)$, with $\hat{\tau^2} = 0.074$ obtained by a hierarchical Bayes model as discussed in (ref). This is reminiscent of the method-of-moments weights in Equation (ref).
Accounting for differences in precision across estimates undoes -- and in fact reverses -- the estimated positive effect of targeting the long-term unemployed. The intuition is most easily seen by looking at the study-level estimates and standard errors for programs targeting this group. Figure (ref) shows the estimated treatment effect and 95% confidence interval for every study targeting the long-term unemployed. A vertical dark line marks the precision-weighted average treatment effect for all studies targeting the reference group (UI recipients). Two studies report treatment effects of at least 30 percentage points -- an enormously large effect relative to much of the literature -- but their estimates are imprecise, with standard errors larger than for many of the other studies. The studies with small or negative effects, by contrast, have small standard errors (or no standard errors at all). A precision-weighted regression gives little weight to the first group and substantial weight to the second; this is why, in Figure (ref), the coefficient for targeting the long-term unemployed becomes negative under WLS. The same qualitative pattern emerges when we additionally account for latent heterogeneity using $w_i = 1/(\tau^2 + \sigma_i^2)$.
Regressions such as equation ((ref)) should be interpreted as predictive models for the effect $\theta_0$ in a new problem instance characterized by covariates $X_0$, with posterior given by Equation (ref). The predicted effect is a weighted average of the previously reported estimates, with weights determined by the prior covariance kernel $K(\cdot,\cdot)$ and by the standard errors $\sigma_i$. This connects directly to the question of external validity: how confident can a policymaker or researcher be that the conclusions of the existing literature carry over to a new setting?
To illustrate, consider two hypothetical new contexts in which a policymaker is contemplating a training program for recipients of mixed gender evaluated at the long-term horizon, again using the card2018works data:
We focus on the long-term unemployed for continuity with the meta-regression discussion above, and compare to UI recipients because they are the most commonly observed recipient group in the meta-analysis sample.
For each context, we compute the posterior mean prediction for $\theta_0$, the posterior standard deviation, and the implicit weights $w_i$ on each prior estimate, under three models that differ only in the prior covariance kernel: (i) weighted least squares with random effects; (ii) a squared-exponential GPP with a long length scale $\ell$, smoothing more aggressively across the covariate space; and (iii) a squared-exponential GPP with a short length scale $\ell$, putting most weight on very close neighbours.
\paragraph{Predictions} Prior hyperparameters are set to sensible defaults based on the data. Covariates are pre-normalised (mean 0, variance 1) before kernels are evaluated. The sample variance of $\hat\theta_i$ is $\widehat{\mathrm{Var}}(\hat\theta) \approx 0.012$ and the median sampling variance is $\widetilde{\sigma}^2 \approx 0.001$, giving a signal variance $V \equiv \widehat{\mathrm{Var}}(\hat\theta) - \widetilde{\sigma}^2 \approx 0.011$. For the linear model, we set $\tau^2 = \rho^2 = V/(q+2)$, where $q=\dim(X_i)$.\footnote{The number of covariates $q$ is 10, because the equation (ref) is estimated using the ten dummy variables shown in Figure (ref).} For the squared exponential kernel, we set $\tau^2 = \rho^2 = V/2$. These choices ensure that the unconditional prior variance of $\theta_i$ matches the observed variance. The squared-exponential length scales are $\ell = \sqrt{q}/2 \approx 1.58$ (short) and $\ell = 2\sqrt{q} \approx 6.32$ (long).
Table (ref) reports posterior means and standard deviations under each method, and Figure (ref) displays the weights graphically for the most heavily-weighted estimates (positive or negative). The per-study estimates, standard errors, distances to each test point, citations, and sample sizes for every row displayed in Figure (ref) are reported in Appendix (ref). Across the three specifications, the posterior means for the employment effects range from about 12 to 15 percentage points for Context A and from about 6 to 10 percentage points for Context B. Context A exceeding Context B is a faint echo of the negative LTU coefficient seen in the regression analysis above.
\paragraph{Policy Recommendations} From predicted effects, it is then a short step to policy recommendations using the results from equation (ref). Suppose that the decisionmaker is willing to pay \$10,000 for each additional employed person and the program costs \$1,000 to administer per-participant, regardless of who is served. The estimates then imply that they would fund the program for UI recipients for whom the expected value is positive ($\$10,000\times0.122-\$1,000>0$) but not for the long term unemployed for whom the expected value is negative ($\$10,000\times0.094-\$1,000<0$) . After receiving such a recommendation, it would be natural for the policymaker (or the researcher generating the recommendations) to want to understand which estimates are most influential for generating the policy recommendation.
\paragraph{Mechanisms} A detailed review of Figure (ref) clarifies the mechanisms underlying the predictions. The left-most column of Figure (ref) shows a reference case that ignores covariates entirely: a fixed-effect inverse-variance-weighted meta-analysis ($\tau^2 = 0$) where each estimate's weight is proportional to $1/\sigma_i^2$. Estimates \#14, \#15, and \#16 receive the most weight because their reported standard errors are roughly an order of magnitude smaller than every other estimate in the sample. Assigning these estimates such a high weight is of course nonsense because they describe programs which are very different from the context of interest.
Focusing first on context A (in sample prediction for UI recipients), Figure (ref) shows that the three models do essentially the same thing, which is to construct an inverse-variance weighted mean of the ten estimates whose covariates exactly match this context. Estimates \#1--\#7 have small standard errors and receive meaningful weight; the other three estimates are sufficiently noisy that they receive little weight.
The differences between the models are most consequential and interesting in context B (out of sample prediction for the long-term unemployed); two patterns are noteworthy. First, meta regression and GP with long $\ell$ place substantial negative weights on a wide set of estimates. With little or no data at the target cell, the linear predictor effectively interpolates from neighbouring cells using contrasts of the form
This assumption is analogous to the parallel trends assumption in a difference-in-difference design.
This re-weighting pattern is visible in Figure (ref) for Context B under meta regression: positive weight on UI training at the long horizon which we used in Context A (Estimates \#1--\#7), positive weight on LTU training at the short horizon (Estimates \#9--\#10), and negative weight on UI training at the short horizon (Estimate \#8). The same mechanism also generates the other negative meta-regression weights visible in the figure.\footnote{For example, negative weight is consistently assigned to Estimate \#16 (Dis/Tr/S, whose unusually small standard error pins down the disadvantaged-inflow indicator and is then subtracted off because the test point is not disadvantaged), To a lesser extent, there are also negative weights on Estimates \#14, and \#15. This occurs for the same underlying reason: because the test point is neither LTU/Other nor Disadvantaged, the model nets out their identifying contributions.} Negative weights aren't a problem in this setting: they are how a linear predictor (or long $\ell$) extrapolates by subtracting off distant test points.
Second, the posterior standard deviations differ dramatically across methods for the out-of-support Context B, even though the point predictions agree. For meta regression, the posterior standard deviation in Context B is essentially the same as in Context A. For the long length scale it grows roughly threefold. For the short length scale it grows by a factor of about seventeen, reflecting that the data are silent about a recipient $\times$ horizon cell with zero observations under a strictly local kernel.
To understand the mechanism by which the posterior standard deviation is so tight in context A (unconditionally) and in context B (when using a long length scale), but is so large in context B when using a short length scale, it is helpful to review Table (ref), which reports the kernel covariances from equation (ref) between the prior studies and the new context. With long $\ell$, the model can confidently predict behavior in the new context. With short $\ell$, the model assumes little covariance between the existing studies and the new context. Instead, the model's prediction is much closer to the precision of the prior before any data is observed (which is the unconditional prior variance equal discussed above under “Predictions”).
How much to extrapolate across the covariate space is ultimately a question about the economic model the researcher is willing to entertain, and the choice of length scale makes that model explicit. A researcher who views employment impacts of the training program as an additive function of the population served and the time horizon will want to borrow heavily across cells. Even if a researcher doesn't view the relationships as literally additive, so long as they are comfortable with extrapolation from far away test points, they can use long $\ell$ GP. Concretely, in the context of this application, these models are suitable if the researcher believes that impacts for UI and LTU differ by one constant, and likewise that long- and short-term impacts differ by a second constant.\footnote{This is reminiscent of the literature on surrogates, where short-term impacts are used to predict long-term outcomes atheyetal2025.} Under such additivity, the long-term effect of training on the long-term unemployed can be reconstructed by adding a population adjustment to the training effect for UI recipients, and a long length scale is appropriate.
A researcher who instead believes that the long-term unemployed respond to training differently --- in either direction --- will prefer a short length scale that declines to make this extrapolation. The structural literature supplies models of both signs. Training may be more effective for the long-term unemployed: in the equilibrium search model of kospentaris2021unobserved, skill loss and reduced search effort account for much of the decline in job-finding among the long-term unemployed, so that human-capital programs do the most for precisely the workers who have lost the most. Or, training may be less effective for them: if employers screen on unemployment duration kroft2013duration, jarosch2019statistical, so that callback rates fall sharply as a spell lengthens, then even retrained long-term unemployed workers may be unable to convert their new skills into offers. Either mechanism implies a training-by-population interaction that the additive predictor cannot represent, and a researcher who finds either one credible should shorten the length scale accordingly.
In Section (ref), we focused on evidence aggregation and extrapolation of findings to new contexts, in the absence of selectivity in the publication process. In practice, however, selectivity is very much a problem for inference based on published evidence: the probability $\bar{d}$ that an empirical finding is published might be a function of that finding itself.
In this chapter, we discuss identification in the presence of selectivity. Such selectivity might be due to decisions by various parties, including researchers (who might engage in p-hacking or specification searching), reviewers and editors (who selectively accept papers, leading to publication bias), as well as due to those conducting a meta-study (who always face the choice of which studies to include). Selectivity does not necessarily reflect improper research methods, but it complicates the interpretation of published findings.\footnote{As discussed formally in whichfindings2018, there are conflicting objectives: validity of conventional inference requires the absence of selectivity. If published findings are to inform subsequent decision-making, however, then surprising findings should be selected, because those are the ones which will have a meaningful impact on decisions. But if there is uncertainty over the validity of some identification strategy or data-source, then surprising findings might be less plausible, and plausibility might require the selection of unsurprising findings.}
Different identification methods differ in terms of their intended purpose and ambition. They might be designed to just test for the presence of selectivity ($\bar{d} \not\equiv 1$). More ambitiously, they might estimate the degree of selectivity and the shape of the selection function $\bar{d}(\cdot)$. They might also provide selection-corrected estimates of average effects $\bar{\theta}$, or more broadly of the distribution $\mu$ of $\theta_i$ across studies. Finally, some methods might provide study-level corrected estimates of the $\theta_i$ themselves.
Different methods for identification furthermore leverage data of varying degrees of richness. First, they might only use the (marginal) distribution of published $Z$-statistics or p-values, and seek to test for jumps or non-monotonicity of the p-curve (density of published p-values). Second, using richer information, they might use the joint distribution of estimates and standard errors (or sample sizes), and seek to leverage a possible dependence between estimates and standard errors to test for the presence of selectivity, or -- more ambitiously -- to identify the degree and shape of selectivity.\footnote{The function $\bar d(\cdot)$ might for instance be U-shaped, if there is a preference for significant results, or inverse U-shaped, if there is a preference for insignificant results, or monotonically increasing/decreasing, if there is a preference for positive/negative effects, etc. See also forking2021.}
Lastly, methods for identifying selectivity might additionally leverage unselected samples of ground-truth estimates. Such estimates might come from research efforts that were pre-registered, for instance by grant agencies, non-academic research labs, or systematic replication studies. Because such non-selected ground truth samples are only rarely available, we focus our review here on the first two approaches, using either just the p-curve, or alternatively the joint distribution of estimates and standard errors.
The results in this section focus on non-parametric identification, and avoid reliance on functional form assumptions. We then discuss both parametric and non-parametric estimation based on these non-parametric identification results in the context of the empirical examples below.
\paragraph{Testing for selectivity: } The derivation in the preceding section assumed that publication was non-selective, so that $D_{i}$ is independent of $Z_{i}$ across studies. In the presence of selectivity, however, we have in general that $\bar{d}(z) = P(D_i=1|Z_i=z)$ is not constant as a function of $z$, and thus, denoting $\tilde{d}(p) = \bar{d}(z)$, $$ g(p) \cdot \frac{\tilde{d}(p)}{E[D_i]} \neq g(p). $$ In words, the density of published p-values deviates from the density $g(p)$ of latent p-values by a factor proportional to the publication probability $\tilde{d}(p)$.
The latent density of $Z$-statistics is given by the convolution $f(z) = \int \varphi(z - \omega) d \nu(\omega)$, where $\nu$ is the distribution of $\omega_{i} = \theta_{i} / \sigma_{i}$. The corresponding density of p-values equals $g(p) = f(z) / \varphi(z)$ for $z= \Phi^{-1}(1-p)$. elliott2022detecting characterize this density, showing that it is smooth and non-increasing. In Appendix (ref) we provide a short proof of these properties.\\
This suggest the following testable hypotheses, which are implied by the null hypothesis of no selectivity:
Tests based of the first of these derived hypotheses (i.e. continuity/smoothness) have power if $\tilde{d}(p)$ is discontinuous at conventional levels. Tests based on the second of these derived hypotheses do not have power if $\tilde{d}(p)$ is also non-increasing: In that case $f(p) = g(p) \cdot \tilde{d}(p)$ is non-increasing, since both $g(p)$ and $\tilde{d}(p)$ are. Typically, researchers prefer smaller p-values, and “significant” results. Whenever that is the case, then $\tilde{d}(p)$ is indeed non-increasing.
There are, however, scenarios where the opposite is true. This includes notably cases where researchers may wish to find an insignificant result, for instance when testing for pre-trends in a difference-in-differences context, and in the context placebo tests for the validity of some research design more broadly.
The approaches that we have just discussed are based solely on the distribution of p-values $P_{i}$ or, equivalently, the distribution of $Z$-statistics $Z_{i} = \hat{\theta}_i / \sigma_i$. These approaches only allow us to test for the presence of selectivity, but not to identify its shape or magnitude, nor to correct published estimates. In meta-studies, however, we typically have not only P-values at our disposal, but instead observe both point estimates $\hat{\theta}_{i}$ and standard errors $\sigma_{i}$. We can therefore leverage their joint distribution for identification, if we are willing to impose some additional assumptions. This also allows us to be more ambitious, going beyond mere tests of the presence of selectivity, to estimate the magnitude and shape of selectivity, and to implement selection corrections.\\
Approaches that leverage this joint distribution typically rely on two assumptions, both of which are substantively restrictive. The first of these is the assumption that estimands and standard errors vary independently across studies, so that
The second is the assumption that selection is only based on p-values,
These two assumption are easily generalized to allow for independence conditional on observed covariates $X_i$; we do so in (ref) below. The first of these two assumptions might be violated if sample sizes in experiments are chosen based on power calculations, and if the researchers preparing these experiments have correct priors regarding how the magnitude of treatment effects varies across studies allcott2015, gechter2024. The second of these assumptions might be violated if selectivity depends on the magnitude of estimates, in addition to p-values, for instance because only certain magnitudes are considered substantively (economically) meaningful by authors or reviewers, regardless of statistical significance.\\
We will discuss three versions of identification based on the joint distribution of $\hat{\theta}_i$ and $\sigma_i$: (1) meta-regressions egger1997bias, (2) the “weighted average of adequately powered” estimates stanley2017finding, and (3) simultaneous estimation of selection function and distribution of estimands publicationbias2019.
All three approaches are based on the same two assumptions: independence and selection based on p-values. The theoretical arguments made below to justify the approaches of egger1997bias and of stanley2017finding do not appear explicitly in the original references, but they rationalize the methods proposed in these papers.
\paragraph{Testing for selectivity: Meta-regressions egger1997bias} Consider the null-hypothesis that there is no selectivity, so that $\bar{d}\equiv 1$. Under this null hypothesis, the independence assumption $\theta_{i} \perp \sigma_{i}|X_{i}$ implies conditional mean independence $E[\hat{\theta}_{i} | \sigma_{i},X_{i}] = E[\theta_{i} | \sigma_{i},X_{i}] = E[\theta_{i} | X_{i}]$. As proposed by egger1997bias, this motivates linear regressions of $\hat{\theta}_{i}$ on $\sigma_{i}$ (and possibly $X_i$, which we omit for notational simplicity in the following), $$ \hat{\theta}_i = \alpha+\beta\cdot \sigma_i + \epsilon _i, $$ or equivalently (this is the version originally proposed by egger1997bias), after dividing both sides by $\sigma_i$,
This latter form of the regression can be motivated as a weighted least-squares estimator that down-weights noisier observations. Note that $\alpha$ and $\beta$ trade places as intercept and slope, after this transformation.
Under the null hypothesis, we obtain $\beta=0$: There should be no systematic relation between the magnitude of estimates $\hat{\theta}_i$ and the size of standard errors $\sigma_{i}$. A conventional $t$-test can be used to test this implication, based on either version of the regression. Rejection of the null that $\beta=0$ can then be interpreted as evidence of selectivity.\\
Note that linearity of the relationship between $\hat{\theta}_i$ and $\sigma_i$ is trivially satisfied under the null, since they are mean-independent. Linearity is however generically violated under the alternative, when the published evidence is selected, so that $\bar{d}\not\equiv 1$: The function
depends nonlinearly on $\tfrac{1}{\sigma}$. This is illustrated by Figure (ref), which plots $E[Z_i|\sigma_i, D_{i} = 1]$, for $D_i = \mathbf 1(Z_i \geq \bar{z})$ with $\bar{z} = 1.96$ and $\theta \sim N(0,1)$, against $1/\sigma$.
This observation is important because many papers aim to interpret $\alpha$ in regression (ref) as an estimate of $\lim_{ \sigma \to 0 } E[\hat{\theta}_{i} | \sigma_{i} = \sigma, D_{i}=1]$ (which in turn is taken to identify the average effect $E[\theta_{i}]$). Such an interpretation is not valid in general. In the example of Figure (ref), any coefficient $\alpha$ between $0$ and $\sqrt{2/\pi}$ could be obtained, depending on the distribution of $\sigma_i$. Linear extrapolation will therefore necessarily lead to biased estimates of the average effect $\bar{\theta}$.
\paragraph{Estimating average effects: Highly powered studies stanley2017finding }
Meta-regressions allow us to construct a test of the null that there is no selectivity (as was the case for our discussion of the distribution of p-values). An alternative goal, in line with traditional meta-studies, is to obtain an estimate of the expectation $\bar{\theta}$ of $\theta_{i}$ across studies, that is, of the average effect. Assume that $\bar{d}(z) \to 1$ as $|z| \to \infty$, which means that findings with large $z$-statistics are always published. By definition, and recalling our assumption that selection depends only on the $Z$-statistic, $$ E[\hat{\theta}_{i} | \sigma_{i} = \sigma, D_{i}=1] = \frac{E[\hat{\theta}_{i} \cdot \bar{d}(\hat{\theta}_{i} / \sigma) | \sigma_{i} = \sigma]}{E[ \bar{d}(\hat{\theta}_{i} / \sigma) | \sigma_{i} = \sigma]}. $$ By assumption, $\bar{d}({\theta} / \sigma) \to 1$ as $\sigma \to 0$, whenever ${\theta} \neq 0$. By the dominated convergence theorem (assuming that $E[|\theta_i|] < \infty$ and $P(\theta_i=0)=0$), this implies $$ \lim_{ \sigma \to 0 } E[\hat{\theta}_{i} | \sigma_{i} = \sigma, D_{i}=1] = E[\theta_{i}]. $$ It follows that the average estimate for “highly powered” studies, with small $\sigma_{i}$, is approximately equal to the average estimand $E[\theta_{i}]$. Put differently, highly powered studies are not subject to selection, because they always yield large $z$-statistics.
This argument justifies the proposed estimator of $E[\theta_{i}]$ in stanley2017finding (see also ioannidis2017power), who suggest to focus on a weighted average of highly powered studies, in order to to estimate $\bar{\theta} = E[\theta_{i}]$. To be empirically feasible, this approach requires the meta-study to include observations in the vicinity of $\sigma_{i}=0$. Put differently, there need to be studies that are highly powered in the sense that $E[\bar{d}(Z_i) | \sigma_i] \approx 1$. In practice, it might well be the case that the latent distribution of $\theta_i$ has mass around $0$, so that even studies with small $\sigma_i$ might not be highly powered. In this case, no reliable estimate of $\bar \theta$ can be formed by using the approach of stanley2017finding. We return to this point in Section (ref).
Note also that the assumption that $\bar{d}(z) \to 1$ as $|z| \to \infty$ requires that both positive and negative estimates are published, whenever $|z|$ is large enough. This excludes one-sided selection, where only positive significant estimates are published, for instance. In the example in Figure (ref) this is not the case. For this example, $\lim_{ \sigma \to 0 } E[\hat{\theta}_{i} | \sigma_{i} = \sigma, D_{i}=1] = 2 \cdot \varphi(0) = \sqrt{2/\pi} \neq 0 =E[\theta_i]$.
\paragraph{Estimating selectivity and effect distributions publicationbias2019} Meta-regressions allow us to test for the presence of selectivity. Highly powered studies allow us to estimate average effects. As it turns out, however, we can achieve considerably more under our assumptions. We can not only test for the presence of selectivity, but we can also identify how much more likely significant estimates are to be published. And we can not only identify the average effect $\bar{\theta}$ when there are highly powered studies, but we can non-parametrically identify the distribution $\mu$ of $\theta_{i}$ across studies -- for instance the variance of effects across studies, or the presence of fat tails. This was proven in publicationbias2019: both the distribution $\mu$ of $\theta$ and the selection function $\bar{d}(z)$ are non-parametrically identified under the assumptions that $\hat \theta_i |\theta_i, \sigma_i^2 \sim N(\theta_i, \sigma_i^2)$, $\theta_{i} \sim \mu$, $\theta_{i} \perp \sigma_{i}$, and $E[D_i =1 | \widehat \theta_i, \sigma_i] = \bar{d}(Z_{i})$.
To provide some intuition for this identification result, we provide both algebra and graphical examples. We first describe how to construct functions of the data which are purged of selectivity, and then describe how to recover the distribution $\mu$ of $\theta$. Algebra and graphical examples are used to illustrate the source of non-parametric identification in the model. Construction of the actual estimators uses maximum likelihood estimation.
Suppose that studies can have only two levels of the standard error: $\sigma_1$ and $\sigma_2$. The distributions of estimates, conditional on either standard error, draw from the same latent distribution of true effects and are affected by the same selection process. They differ only in the fact that $\sigma_2$ is noisier than $\sigma_1$. Assume further (again just for illustration) that the latent distribution $\mu$ is normal.\footnote{publicationbias2019 prove identification in a setting where the distributions of true effects and standard errors are left unrestricted, and selectivity can based on an arbitrary function of the $Z$-statistic.}
Figure (ref) illustrates several distributions in this simplified environment. Two latent distributions are shown in black: in the left panel the mean is zero and in the right panel the mean is three. Both latent distributions have variance, $\tau^2 = 1$. In addition, the figure shows the distribution of estimates drawn from the latent distribution with noise $\sigma = 1$ (published and unpublished) and the distribution of estimates with noise $\sigma = 2$ (published and unpublished).
\captionsetup{font=small} \captionsetup[sub]{font=small,skip=0.0cm}
Consider the density of the $Z$-statistic across studies (both published and unpublished), conditional on the standard error: $$ f(z|\sigma) = \int \varphi\left( z - \tfrac{\theta}{\sigma} \right) d \mu(\theta) $$ Figure (ref) shows the distributions of the $Z$-statistics for the illustrative example. Larger standard errors are associated with smaller $Z$-statistics so the distribution for $\sigma = 1$ is more dispersed than the distribution for $\sigma = 2$.
The density of the $Z$-statistic conditional on the standard error and conditional on observability ($D_{i} = 1$) is then given by $$ f(z|\sigma, D=1) = \frac{\bar{d}(z)}{E[\bar{d}(Z)|\sigma]}\cdot f(z|\sigma). $$
Figure (ref) shows the distributions of the $Z$-statistics in the presence of selectivity. We assume that selectivity can be characterized by two values: the publication rate for significant estimates ($Z > 1.96$), which is normalized to 1, and the publication rate for insignificant estimates ($Z < 1.96$), which is assumed in this example to be 0.5. Both distributions therefore show a discontinuity with a doubling of the density at $Z = 1.96$.\footnote{This plot makes apparent the connection between the meta-studies approach method and the approach based on discontinuities in the distribution of p-values from Section (ref). The latter effectively zooms in on the density close to the discontinuity. The meta-studies method uses all of the information in the density but requires the independence assumption; it also yields much richer results.}
The key identification idea behind the meta-studies model in publicationbias2019 is that if we take the ratio of the density $f(z|\sigma, D=1)$ across these two values of $\sigma$, then the selection probability $\bar{d}(z)$ appears in both the numerator and the denominator so it drops out:
where $const.$ is a constant that does not depend on $z$.
Figure (ref) shows the ratio of the two densities. Notice the ratio of densities in Figure (ref) is the same as the ratio of densities in Figure (ref). Although the densities in Figure (ref) are distorted by selectivity, the ratio of the densities is not distorted. This is because the density function for $\sigma = 1$ jumps by a factor of $\frac{1}{0.5} = 2$ at the cutoff, and the density function for $\sigma = 2$ also jumps by a factor of two at the cutoff. Having established how to construct functions of the data which are purged of selectivity, we next turn to how parameters of the distribution can be identified using these functions of the data.
The relative distribution of $Z$-statistics---encapsulated by the density ratio---provides a “signature” which can be uses to identify the mean of the latent distribution. In Figure (ref), when the latent mean is zero, the density ratio of the $\sigma=1$ distribution to the $\sigma=2$ distribution is greater than 1 when $z<-1.4$ or $z>1.4$. In contrast, when the latent mean is three, the density ratio of the $\sigma=1$ distribution to the $\sigma=2$ distribution is less than $1$ for a much larger region (from $z>-4.7$ to $z<2.7$). Put simply, when the latent mean is positive, estimates from the precise distribution will more often be associated with large positive $Z$-statistics than estimates from the imprecise distribution. The larger the region where the density ratio is less than 1, the greater is the latent mean. This is the “reduced form” signature of a positive latent mean in the presence of selectivity. We show this same signature exists in data from cohen2025disemployment in (ref).
The result of publicationbias2019 formalizes this intuition, showing the distribution $\mu$ of $\theta$ is uniquely pinned down by the ratio shown in (ref). Once $\mu$ is identified, we can recover $\bar{d}$ (up to a multiplicative constant) from $$ \bar{d}(z) = const. \cdot\frac{f(z|\sigma, D=1)}{f(z|\sigma)}. $$
Since identification of both $\mu$ and $\bar{d}$ is non-parametric, they could in principle be estimated non-parametrically (e.g. using GMM, as discussed in the appendix of publicationbias2019). For moderate sample sizes (number of estimates in the meta-study), however, parametric estimation tends to be better behaved.
We investigate the extent of selectivity across all of our running applications. The nature of selectivity depends upon researchers' priors and preferences about the size of the causal effect of interest. In all four of our applications, economic theory or nudge design predicts a positive effect of the treatment (e.g. higher re-employment, longer unemployment, higher consumption, increase in outcome targeted by the nudge). The concern regarding selectivity is therefore that researchers and the publication process will select for positive significant effects (rather than positive and negative significant effects).
We report results from three types of tests for selectivity discussed above. For elliott2022detecting, we test whether there is a discontinuity in the density of one-sided p-values and $Z$-statistics (a monotonic transformation of p-values) at $p = 0.025$ (corresponds to $Z = 1.96$) using a local polynomial density estimator cattaneo2020simple. Estimates are shown in columns 4 and 5 of (ref).\footnote{We show histograms of p-values and discuss sensitivity of our results to bandwidth selection in Appendix (ref).} Second, for egger1997bias, we estimate equation (ref) and report $\beta$, where $\beta\neq0$ indicates the presence of selectivity. Third, we implement estimators using the model in publicationbias2019, both assuming a latent $t$ distribution and making no assumption about the shape of the latent distribution. Given the discussion above about selectivity being specifically likely to generate positive significant effects, we report estimates for $P(\text{publish} | Z < 1.96) / P(\text{publish} | Z > 1.96)$.
We find widespread evidence of selectivity. We find evidence of selectivity---in the sense that we can reject the null hypothesis of no selectivity---for three out of four applions using the tests from publicationbias2019 and the test fromegger1997bias.\footnote{card2018works reports a finding of no selectivity based on the test in egger1997bias; we replicate that finding,} The evidence using the density discontinuity approach from elliott2022detecting is nuanced. In every case, the estimate implies that there is substantial selectivity at the cutoff but we are only once able to reject the null that the density ratio at the cutoff is 1. The fact that we both estimate substantial selectivity and are unable to reject the null indicates the test for discontinuity is low-powered in our contexts.\footnote{elliott_power_2025 show that in some contexts the power of discontinuity tests can be quite low. So a failure to reject the null of no selectivity should be interpreted as dispositive. We also find that testing for $f(p)$ increasing in $p$ is low-powered (another concern raised by elliott_power_2025 ).} Using the MLE method with a latent $t$ distribution, we find that the probability of reporting insignificant results relative to significant ones varies from 5% to 28% depending on the application.
A finding of selectivity may be surprising in the context of RCTs (nudges and cash transfers) because trial pre-registration should in principle limit the selective reporting of findings. One possibility is that selectivity emerges not in terms of which studies report results, but which outcomes are reported or emphasized by study authors. dellavigna2022rctstoscale find evidence of selective reporting in terms of the most significant $Z$-statistic (as opposed to all $Z$-statistics). Although crosta2024unconditional works to standardize outcomes as much as possible, there is still a range of consumption definitions which researchers can report. Sometimes food is excluded and sometimes it is included, sometimes durables are excluded and sometimes they are included, sometimes the time horizon is one month and sometimes it is one year. Another possibility is that the identifying assumption of orthogonal estimates and standard errors is violated through power calculations, where there are some studies which are projected to have bigger effects for that subpopulation, but estimates for that subpopulation rely on smaller trials.\footnote{As one concrete example of this, many of the studies with the highest standard errors have framing which is tied to specific life circumstances (e.g. small business entrepreneurship or children). One might plausibly think that there are larger effects of cash in this subpopulation (indeed this is the motivation for targeting in the first place). }
We next ask whether selectivity affects the distribution of the latent estimates across our applications in Table (ref).
We use two approaches to correct for selectivity.\footnote{We do not use the correction approach from equation ((ref)) because, as shown in Section (ref), it is not identified; we do not use the correction approach from MLE nonparametric since we do not believe it is robust in our settings.} Our first approach analyzes highly powered studies using the approach proposed in ioannidis2017power. We estimate $\hat{\alpha}$ from equation ((ref)) and then take the subset of estimates with $\sigma_i<|\hat{\alpha}|/2.8$. ioannidis2017power say that this approach filters to the studies which will capture the true latent mean (of $\hat{\alpha}$) with 80% power. They then calculate a precision-weighted mean among this subset of “highly-powered” studies. The highly-powered studies method (at least as implemented following the procedure in ioannidis2017power) only delivers a meaningful subset of studies to analyze in two out of four meta-analyses. The fact that there are very few studies which meet this criteria motivates alternative selection correction methods that use all the available data. In our second approach, we use the model from publicationbias2019.
Correcting for selectivity using this method reduces the mean of the latent distribution $\bar\theta$ in our four applications. More specifically, the ratio of the corrected mean to the simple mean ranges from 12% to 21% depending on the application.\footnote{Our estimates imply that the latent distribution for dellavigna2022rctstoscale is Cauchy. The mean is undefined for a Cauchy distribution. In this case, $\bar\theta$ can be interpreted as the median.} In the context of selection oriented towards positive significant findings (which we believe is a good description of the four meta-analyses that we study), this is precisely what one would expect to find. The distribution of published findings---which features many positive significant findings but has “deleted” some of the positive insignificant and negative insignifcant findings---will have too high of a mean. Selectivity can therefore lead an unsophisticated literature review to overstate the average causal effects of the policies studied.
Table (ref) also reveals that the precision weighting and incorporating fat tails from Section (ref) provide a partial back door approach to selection correction. If the true data-generating process involves selectivity in favor of positive significant findings, then estimates with large standard errors will on average be more biased than estimates with small standard errors. Precision-weighting reduces the influence of these more-biased estimates. Our findings also reflect this pattern: the precision-weighted means (shown earlier in Table (ref) and repeated for convenience in Table (ref)) are smaller than the simple means in all four applications.\footnote{This has the interesting implication that selection correction will not always change the conclusion of a meta-analysis. As one concrete example, crosta2024unconditional focus on MLE estimates which precision weight each estimate (but do not correct for selectivity). In that case, the simple mean for the estimated average effect of \$1 on monthly consumption is \$0.12, the estimated effect when using MLE without correcting for selectivity is \$0.03, and the estimated effect of using the MLE with correcting for selectivity of \$0.02. So even without correcting for selectivity, precision weighting generates similar conclusions about the mean of the effect distribution to a selection correction.} This is what we would expect to see in the presence of selectivity; in the absence of selectivity, the precision-weighted means could plausibly be larger or smaller than the simple means.
Last but not least, it is useful to note that the empirical estimates from the MLE latent $t$ model falsify an identifying assumption from the highly-powered approach. The problem, which we describe in the theory section, is that if there are some latent effects $\theta_i$ which are zero, then no standard error will be sufficiently small as to qualify as “highly powered.” Across all applications, the distributions implied by the estimates for $\bar\theta$ and $\tau$ imply that there exists a meaningful distribution of latent effects in the neighborhood of zero.
Many meta-analyses want to combine covariates and selectivity -- that is, they want a predictive model for new problem instances characterized by covariates that also accounts for selectivity.
Our discussion in (ref) sidestepped study level covariates. Covariates can, however, be incorporated either by partitioning the observed studies and then separately implementing any of the selectivity methods\footnote{See for example brodeur2020methods, who consider p-curves separately depending on the methods employed in different empirical studies.} or by directly building covariates into the models used. The former approach requires sample sizes larger than those present in our empirical examples so we will focus on the latter.
\paragraph{Likelihood} Recall the assumptions stated in (ref), where we assumed a sampling distribution of $\hat \theta_i |\theta_i, \sigma_i^2,X_{i} \sim N(\theta_i, \sigma_i^2),$ and selection based solely on $Z$-statistics $Z_i = \hat \theta_i/\sigma_i$ according to $P(D_i =1 | \hat \theta_i, \sigma_i, X_i) = \bar{d}(Z_{i}).$
Suppose now additionally, as in the meta-regression model (equation (ref)), that $$ {\theta}_i| X_i =x \sim N(x\cdot \beta, \tau^2). $$ Then, again as before, $\hat{\theta}_{i} |X_{i}, \sigma_{i} \sim N(x \cdot \beta, \tau^{2} + \sigma_{i}^{2}),$ and thus $$Z_{i} |X_{i}, \sigma_{i} \sim N(x \cdot \beta / \sigma_i, \tau^{2}/\sigma_{i}^{2} + 1).$$ We get the conditional likelihood of published $Z$-statistics by Bayes rule (this implicitly also conditions on model parameters $\beta,\tau^{2}$ and $\bar{d}$), $$ f(Z_{i} |X_{i}, \sigma_{i}, D_i = 1) = \underbrace{ \frac{1}{\sqrt{\frac{\tau^{2}}{\sigma_i^{2}}+1 }} \varphi\left(\frac{Z_{i}-X_{i}\cdot\beta}{\sqrt{\frac{\tau^{2}}{\sigma_i^{2}}+1 }} \right) }_{ f(Z_i|X_i, \sigma_i) }\cdot \underbrace{ \frac{\bar{d}(Z_{i})}{E[\bar{d}(Z)| X_{i}, \sigma_{i}]} }_{ \frac{P(D_i=1 |Z_i, \sigma_i,X_i)}{P(D_i=1|X_i, \sigma_i)} }. $$ If we consider a parametric specification for $\bar{d}$ of the form $\bar d(z) = \mathbf{1}(z>1.96) + \gamma\cdot\mathbf{1}(z \leq 1.96)$ (selection for positive significant effects), then the denominator in the last fraction equals $$E[\bar{d}(Z)| X_{i}, \sigma_{i}] = 1 + (\gamma-1) \cdot\Phi\left(\frac{1.96-X_{i}\cdot\beta}{\sqrt{\tfrac{\tau^{2}}{\sigma_i^{2}}+1 }} \right). $$
\paragraph{Prior and posterior} We can combine this likelihood with a prior for the hyper-parameters $\beta,\tau^{2}$ and $\gamma$ to get a complete hierarchical Bayes specification. We can then use Hamiltonian Monte Carlo (as implemented in Stan) to sample from the posterior. Using the posterior mean $\bar{\beta}$ of $\beta$, we can furthermore make predictions for new instances, via $X_{0}\cdot\bar{\beta}$.
For our empirical application, we assume a normal prior for $\beta$, $\beta \sim N(0, \Sigma)$, as in the Ridge regression models discussed earlier, a half-normal prior for $\tau^2$ with large variance, and a lognormal prior for $\gamma$, $\log(\gamma) \sim N(0,1)$.
The meta-analyses that want to combine covariates and selectivity include three out of four of our applications: cohen2025disemployment, dellavigna2022rctstoscale, card2018works.\footnote{The first article studies the combination using Bayesian Model Averaging, the second article studies the combination by using a complex simulation procedure (Table VI), and the third article is a priori interested in both questions, but did not find evidence of selectivity and therefore views covariate-only models as sufficient.} cohen2025disemployment in particular is interested in the optimal replacement rate for unemployment insurance. It finds that higher baseline replacement rates ($X_i$) are associated with a higher elasticity of unemployment duration with respect to benefits ($\hat\theta_i$). In the model they study, solving for the optimal replacement rate requires predicting how the elasticity varies with baseline replacement rate.
(ref) builds up to the model described in (ref). We start with a multivariate OLS regression. cohen2025disemployment find that estimates for the elasticity of unemployment duration with respect to replacement rate are larger in contexts with a larger baseline replacement rate. We replicate this finding by running the minimal regression needed to investigate this interaction effect:
where “RR” is replacement rate and “PBD” is potential benefit duration.
Then, we account for differing precision by switching to weighted least squares with weights $w_i = \frac{1}{\sigma^2}$. This leads to large changes in the coefficients and in particular shows no difference in elasticity with respect to RR for changes in the baseline RR. Accounting further for heterogeneity through weighting by $w_i = \frac{1}{\sigma_i^2 + \tau^2}$ (with the $\tau^2$ obtained by estimating the hierarchical Bayes model without a selection correction) leads to more equal weights across studies; generally this will lead to coefficent estimates that fall between those of OLS and WLS with $w_i = \frac{1}{\sigma^2}$.\footnote{In a multivariate regression this is not a rule -- as demonstrated by this relationship failing with covariate “PBD x baseline RR” -- because changes in the values of other coefficients will change the coefficients of collinear covariates. The hierarchical Bayes model without a selection correction almost exactly replicates this WLS with $w_i = \frac{1}{\sigma_i^2 + \tau^2}$.} We get identical results by running the hierarchical Bayes model without a selection correction.
The final model accounts for selection. We find that when the replacement rate increases by 10%, the elasticity increases by 8.6%. The hierarchical Bayes model without a selection correction would estimate an elasticity increase of 5.6%. The selection-corrected model implies an optimal replacement rate of 32% in the US context studied by cohen2025disemployment if all other factors except the slope of the elasticity vs. baseline replacement rate are held constant. The uncorrected model implies an optimal replacement rate of 29%.
This section collects our recommendations into a single pipeline. The pipeline runs in six stages, summarized in (ref): (0) scoping and data, (1) aggregation, (2) prediction with covariates, (3) selectivity, (4) the combined model, and (5) decisions.
\paragraph{Stage 0 -- Scope, collect, assess credibility.} Fix the estimand $\theta_0$, the search procedure, and the coding rules and write down the proposed procedure. This matters because apparent selectivity can be an artifact of which outcome was coded rather than which study was published, as in the cash-transfer application ((ref)). AI tools make extraction cheap; see cook2026aimeta for pitfalls. The methods in this paper make sense when there are at least three prior estimates, so you should stop if you cannot find three prior estimates.
At this stage, also decide whether to collect covariates and which covariates to collect, since it determines the model in Stage 3. Do not discard lower-quality studies. Instead, collect measures of study quality as part of the covariates (and return to this question in step 3). This could be something relatively objective like the research design (e.g. DiD, RD, RCT), a subjective assessment of the credibility of the paper's identifying assumptions, or both.\footnote{Although we recommend this practice, we did not implement it as part of our empirical examples in this article. That is because the notion of study quality is difficult to treat in a consistent way across very different applications. We advocate for using the meta-analyst's assessment of quality (rather than say journal quality, which could be related to selectivity forking2021.}
\paragraph{Stage 1 -- Aggregate without covariates or selection.}
Ask whether the dispersion of the $\hat\theta_i$ exceeds what their standard errors imply ((ref)). Report a precision-weighted mean ((ref)), not a simple mean, and form shrunken study-level estimates ((ref)). Whether you proceed to the subsequent steps depends on the number of estimates you have collected. The subsequent steps need a soft floor of roughly thirty estimates (irsova2024meta).\footnote{We take the upper bound from irsova2024meta who suggest at least 10 primary studies with 30 total estimates.}
Next, if you have enough studies, fit and compare three latent distributions ((ref)): normal , $t$ for heavy tails, and nonparametric via nonparametric maximum with Tweedie posterior means ((ref)). Read off whether tails are heavy (small $\nu$) and whether $\mu$ is skewed -- these choices materially move the estimated mean and shrink $\tau^2$.
\paragraph{Stage 2 -- Predict with covariates.} Standardize covariates by their empirical standard deviation. For a context inside the support of the covariates, use a meta-regression estimated as random effects WLS ((ref)), equivalently the Gaussian-process posterior mean ((ref)); the equivalence makes explicit that the prediction is a transparent weighted average of the $\hat\theta_i$. Although we emphasize the importance of distributions that allow for fat tails in stage 1, this becomes less important once covariates are available; this is because it is likely that the distribution is assumed to be normal after conditioning on the prediction $X_i\beta$ or $\bar{\theta}(X_i)$, where the covariates $x$ can explain why a specific estimate $\hat\theta_i$ is an outlier.
For extrapolation outside the support, use a Gaussian process prior with a squared-exponential kernel (equation (ref)). Make an active decision about whether to use a long or short length scale. The short length scale will lead to a much larger posterior standard deviation. This choice should be informed by your economic model of how outcomes are determined (see discussion in Section (ref) as an example).
A second use of covariates is to incorporate study quality into the prediction problem. If you are interested in the likely effects of a policy, do the prediction exercise where $\theta_0$ is defined as the highest quality strata of study.\footnote{As one example, cohen2025disemployment report predictions for a regression discontinuity design when doing policy analysis.} The researcher's choice of prior will decide how much consideration to give to the studies from the lower-quality strata. Another option is to do the pipeline for just the high-quality studies; this is implicitly what is already done in the context of the small-sample meta-analysis (Tables (ref), (ref), and (ref)).
For any such prediction model for $\theta_i$ given $X_i$, it is important to keep in mind that (i) the resulting coefficients are not causal and (ii) any coefficients in a multivariate regression capture variation holding all other regressors constant.
\paragraph{Stage 3 -- Test and correct for selectivity.} Assess whether you believe that equations ((ref)) and ((ref)) are satisfied.\footnote{As an example of where they may not be satisfied, dube2024minwage find in the context of minimum wage own-wage elasticities, “credible designs may use less minimum wage variation, resulting in lower precision, but they may also have less bias”. In randomized control trials, researchers have priors about the size of the causal effect they use in power calculations to choose sample sizes, leading to larger expected causal effects being systematically less precise. See chen2025precision for an example of a model with co-dependence (but is not about meta-analysis)} If they are not satisfied, then you are limited to the a $p$-curve density-discontinuity test elliott2022detecting. However, it is important to know that discontinuity tests can have low power elliott_power_2025, so non-rejection may not imply the absence of selectivity. If equations ((ref)) and ((ref)) are satisfied, you can also use the egger1997bias meta-regression, and the relative publication probability from publicationbias2019.
For correction, we emphasize that two commonly used methods are invalid. The Egger intercept is a valid test but not a valid bias correction -- under selection, $E[Z_i\mid\sigma_i]$ is nonlinear in $1/\sigma_i$ and linear extrapolation to $\sigma_i=0$ can return almost any value ((ref)). The highly-powered-studies estimator stanley2017finding, ioannidis2017power is valid only when such studies exist; when $\mu$ has mass near zero they may not, as in two of our four applications ((ref)). Our workhorse is the publicationbias2019 model, which jointly estimates $\mu$ and the selection function from the joint distribution of $\hat\theta_i$ and $\sigma_i$ ((ref)).
One note of caution regarding the correction method of publicationbias2019 is that it can deliver results which are fragile to assumptions about the distribution of latent effects as well as the chosen thresholds for selective reporting. Although we report only one selectivity-corrected MLE model here for brevity, we recommend that authors assess sensitivity to both assumptions about the distribution of latent effects as well as chosen thresholds. See cohen2025disemployment for an example of such sensitivity tests.
\paragraph{Stage 4 -- Combine covariates and selectivity.} Build covariates into the selection model: the linear meta-regression $\theta_i\mid X_i\sim N(X_i\cdot\beta,\tau^2)$ combined with a step-function $\bar d$, estimated by Hamiltonian Monte Carlo in Stan. We recommend doing this with the \href{https://github.com/wwiecek/baggr/tree/master}{baggr} package. Build it up in visible steps -- OLS, WLS, random effects WLS, hierarchical Bayes without and then with selection ((ref)) -- so the reader sees what each ingredient does.
\paragraph{Optional Stage 5 -- From estimates to decisions.} For a binary choice, carry the prediction through to the rule rather than stopping at a coefficient. Under welfare $U(A_i,\theta_i)=A_i\theta_i$, choose $A_i=\mathbf 1\!\left(E[\theta_i\mid \hat\theta_i,\sigma_i,X_i]\ge 0\right)$ for an evaluated site and $A_0=\mathbf 1\!\left(E[\theta_0\mid X_0,\hat\theta_{1:n}]\ge 0\right)$ for a new context (equation (ref)).