Sensitivity and identification quantification by a relative latent model complexity perturbation in Bayesian meta-analysis.

In recent years, Bayesian meta-analysis expressed by a normal-normal hierarchical model (NNHM) has been widely used for combining evidence from multiple studies. Data provided for the NNHM are frequently based on a small number of studies and on uncertain within-study standard deviation values. Despite the widespread use of Bayesian NNHM, it has always been unclear to what extent the posterior inference is impacted by the heterogeneity prior (sensitivity S ) and by the uncertainty in the within-study standard deviation values (identification I ). Thus, to answer this question, we developed a unified method to simultaneously quantify both sensitivity and identification ( S - I ) for all model parameters in a Bayesian NNHM, based on derivatives of the Bhattacharyya coefficient with respect to relative latent model complexity (RLMC) perturbations. Three case studies exemplify the applicability of the method proposed: historical data for a conventional therapy, data from which one large study is first included and then excluded, and two subgroup meta-analyses specified by their randomization status. We analyzed six scenarios, crossing three RLMC targets with two heterogeneity priors (half-normal, half-Cauchy). The results show that S - I explicitly reveals which parameters are affected by the heterogeneity prior and by the uncertainty in the within-study standard deviation values. In addition, we compare the impact of both heterogeneity priors and quantify how S - I values are affected by omitting one large study and by the randomization status. Finally, the range of applicability of S - I is extended to Bayesian NtHM. A dedicated R package facilitates automatic S - I quantification in applied Bayesian meta-analyses.

INTRODUCTION

To establish reliable evidence regarding the benefits of various healthcare interventions, we need statistical methods that can cope with varying results from different studies and allow for a combined evaluation of evidence. Under these circumstances, meta‐analysis is the statistical method of choice in evidence‐based medicine, because it shifts the focus away from the evidence in a single study toward a combined view of evidence provided by several different studies. Such analysis takes into account the between study heterogeneity that, when ignored, can result in overoptimistically small standard errors of estimates, meaning that the significance of an estimated effect may actually be invalid. Frequently, however, only a very small number of studies are available for meta‐analysis, and thus, in these cases, the results can also be misleading due to unobserved heterogeneity. The Bayesian methodology mitigates this problem by incorporating priors on heterogeneity standard deviation into a normal–normal hierarchical model (NNHM).

However, the heterogeneity in a meta‐analysis based on studies can only be well understood when both the between‐study heterogeneity and the within‐study standard deviation values are simultaneously considered (Hardy & Thompson, 1998). The mutual interplay between the between‐study heterogeneity and the within‐study standard deviation values is also active in the Bayesian setting and impacts the posterior estimates (Gustafson et al., 2006). However, this interplay is prone to two sources of uncertainty. While it is well known that the choice of the prior on the heterogeneity parameter is uncertain, it seems to be less well known that the within‐study standard deviation values are actually uncertain estimates. Both sources of uncertainty act simultaneously and can affect posteriors. Therefore, we need a method that can simultaneously quantify the impact of both sources of uncertainty.

Although only two studies are necessary for a Bayesian NNHM to be applicable, in cases when only a small number of studies are supplied for a meta‐analysis, the impact of the heterogeneity prior on the posterior can become a concern. In many cases, various selection criteria and exclusion of old and outlying studies can generate subgroups of the original data set, thus further reducing the number of studies included (Bjordal et al., 2004; Crins et al., 2014; Smith et al., 1995). Indeed, over of meta‐analyses extracted from the Cochrane Database of Systematic Reviews (CDSR) are based on only two studies and nearly 75 are based on five or fewer studies (Davey et al., 2011; IntHout et al., 2015; Turner et al., 2012). When only such a small number of studies are included in a Bayesian NNHM, the sensitivity of the posterior to the heterogeneity prior can become an issue.

Bayesian meta‐analysis expressed by an NNHM is also used to analyze studies when the data are too limited to provide a precise estimate of the within‐study standard deviation . There are studies that contain only two participants and at least 25 of the studies of the CDSR contain less than 50 participants (Davey et al., 2011; Turner et al., 2012). When studies are not sufficiently large, then the estimation of the within‐study standard deviation is unstable (Hoaglin, 2015). Therefore, when such studies with only a small number of participants are included in a Bayesian NNHM, the sensitivity of the posterior to the uncertainty in the within‐study standard deviation values is a concern.

Practical applications of the Bayesian NNHM thus need to be aware of both the sensitivity of the posterior to the heterogeneity prior (called here sensitivity ) and the sensitivity of the posterior to the within‐study standard deviation uncertainty (called here identification ). Moreover, these concerns are intertwined and cannot be separately assessed, thus necessitating a simultaneous two‐dimensional sensitivity‐identification quantification. However, we currently still lack a systematic account of how sensitivity and identification affect the posterior inference in a Bayesian NNHM.

As it is impossible to elicit a unique prior (Gelman et al., 2017), a method is needed that perturbs the specific prior and assesses the impact of prior perturbation on posterior inference. The research area for formal sensitivity quantification is very broad and comprises, for example, the robust Bayesian analysis (Berger & Berliner, 1986; Berger et al., 2000; Ríos Insua et al., 2000), local robustness (Gustafson, 2000) and the epsilon‐local sensitivity (Roos et al., 2015). The local sensitivity approach deals with the rate of change of posteriors with respect to prior perturbation, usually using differential calculus to quantify local sensitivity.

Identification is a concern in a Bayesian setting, because a model expressed by a nonidentified likelihood can be still estimable within the Bayesian framework, given prior assumptions. This means that in the Bayesian setting, the posterior inference for nonidentified parameters is not impacted by the data but is completely determined by the priors. This can be problematic for scientific inference and policy decisions, which endeavor to be determined by data, because in a Bayesian hierarchical model, it is difficult to know which parameters are informed by the data and which are nonidentified.

The broad and complex research area of formal identification investigation includes Bayesian identifiability (Gelfand & Sahu, 1999; Xie & Carlin, 2006) and Bayesian inference for partially identified models (Gustafson, 2015). Identification addresses the impact of the number of observations in a data set (likelihood) on the posterior inference. Therefore, some approaches to identification explicitly change (increase or reduce) the number of observations in the analysis. For example, in the context of classical statistical analysis, both cloning and subsetting of observations have been suggested (Jiang, 2013; Lele et al., 2010; Normand, 1999). Neither data cloning nor the subset argument is, however, applicable if there are only two studies included for a medical meta‐analysis. Thus, in order to quantify identification in a medical Bayesian meta‐analysis, the likelihood must be perturbed such that its impact on the posterior system is changed, but without explicitly adding or removing observations.

To address these concerns, we develop a principled approach to combined sensitivity and identification quantification in a Bayesian NNHM. The novel two‐dimensional sensitivity‐identification (‐) measure addresses sensitivity and identification in a Bayesian NNHM for all model parameters on an equal footing. We keep the prior for the main effect fixed and focus on perturbations of the heterogeneity prior for and the within‐study standard deviation . Using the approach to aligning heterogeneity priors for a Bayesian NNHM proposed by Ott et al. (2021) as a starting point, we suggest systematic likelihood and heterogeneity prior perturbations that are based on perturbations of the relative latent model complexity (RLMC). While sensitivity quantifies the impact of the heterogeneity prior controlled by its scale parameter, identification is concerned with the impact of the data (‐uncertainty) on the marginal posterior inference. We provide ‐ estimates for one light‐tailed half‐normal (HN) and one heavy‐tailed half‐Cauchy (HC) heterogeneity prior and three RLMC targets. Finally, we provide the R code for automatic ‐ quantification in an R package si4bayesmeta on R‐Forge.

We apply the proposed framework for simultaneous ‐ quantification to three medical case studies: historical data for a conventional therapy (Table 1), data with one large study included and excluded (Table 2), and two subgroup meta‐analyses specified by their randomization status (Tables S4 and S5 of the Supporting Information). Providing concrete ‐ estimates helps to pinpoint the severity of the impact of both the heterogeneity prior and the uncertainty on posterior estimates, which allows the impact of HN and HC heterogeneity priors to be assessed.

TABLE 1

Data for the mean pocket reduction (MPR) after a conventional nonsurgical periodontal treatment (Zaugg et al., 2014): within‐study sample size (total), MPR (mean), standard deviation (sd), and data supplied for the Bayesian meta‐analysis: MPR () and (). The reference standard deviation for this data set is (Equation (4))

	Conventional			y	σ
Study	Total	Mean	sd	MPR	SE(MPR)
1	10	0.28	0.32	0.28	0.10
2	11	2.25	2.47	2.25	0.74
3	14	3.00	1.57	3.00	0.42
4	10	0.80	0.86	0.80	0.27
5	42	0.47	0.11	0.47	0.02
6	84	1.20	1.33	1.20	0.15
7	17	0.10	0.30	0.10	0.07
8	24	1.62	0.58	1.62	0.12
9	8	1.18	0.81	1.18	0.29
10	30	0.84	0.55	0.84	0.10
11	30	1.67	2.61	1.67	0.48
12	19	1.95	1.94	1.95	0.45
13	27	0.90	2.00	0.90	0.38

TABLE 2

Data for the association between diabetes mellitus (DM24) and the severe course of Coronavirus Disease 2019 in studies (Kumar et al., 2020): total within‐study sample size (), (), and (). The reference standard deviation of the DM24 data set is (Equation (4)). All studies with exception of the largest study based on registry data from 6637 patients (Study 20) had good quality. Study 20 was only classified as a fair quality study. Removal of Study 20 leads to a smaller data set (DM23) based on studies with the reference standard deviation (Equation (4))

Study	n Total	y log(OR)	σ SE(log(OR))
1	393	0.200	0.244
2	140	0.261	0.519
3	116	0.305	0.516
4	221	0.384	0.487
5	214	0.419	0.395
6	476	0.463	0.320
7	161	0.588	0.862
8	273	0.641	0.505
9	548	0.651	0.245
10	597	0.798	0.398
11	124	0.925	0.538
12	487	1.148	0.463
13	1099	1.157	0.250
14	1012	1.211	0.453
15	155	1.301	0.667
16	298	1.303	0.499
17	201	1.471	0.503
18	119	1.477	0.648
19	138	1.520	0.581
20	6637	1.528	0.109
21	135	2.186	0.698
22	167	2.315	0.666
23	123	2.333	0.778
24	120	4.056	1.479

In the second section of this paper, we review two challenging and motivating case studies. In Section 3, we introduce the general theoretical framework supporting our method, specify the numerical computation of both raw ‐ and relative ‐ estimates, and introduce the R package si4bayesmeta. Results presented in the tables in Section 4 demonstrate the relevance of both ‐ and ‐ estimates for medical case studies. For the sake of brevity, some relevant material is then included only in the Supporting Information. Primarily, the unified ‐ approach is extended to cover Bayesian normal‐t hierarchical model (NtHM) fitted by MCMC in JAGS. The Bayesian NtHM addresses the important case when the assumption of the normality of random effects is not satisfied and thus provides a robustification of the popular but simple Bayesian NNHM against outliers in the data. This extension proves useful in one medical case study and in a simulation study. Moreover, we demonstrate how the randomization status of two subgroup meta‐analyses affects the ‐ estimates. We conclude with a discussion in Section 5.

DATA

This section reviews two challenging and motivating case studies. The historical data set for a conventional periodontal therapy is focused on the mean pocket reduction (MPR), which is a continuous outcome. The second application, which has a binary outcome, deals with the association between diabetes mellitus and the severe course of Coronavirus Disease 2019 (COVID‐19). It consists of two data sets, one of which is a subset of another data set from which one large study has been excluded. One additional case study, which consists of two subgroup meta‐analyses specified by their randomization status is provided in Section 3 of the Supporting Information.

Continuous outcome: Mean pocket reduction

Periodontitis is a complex and increasingly prevalent infective disease that affects the tissue surrounding teeth (GBD, 2016). As a clinical consequence, tissue inflammation leads to pocket formation, bone resorption, tooth loosening, and tooth loss. The conventional nonsurgical periodontal therapy involves removal of dental plaque and calculus (scaling), followed by smoothing (planning) of the affected surfaces of the roots. Usually, the success of this therapy is quantified by the reduction of the pocket depth after treatment. High values of MPR indicate benefit to the patient. Table 1 shows historical, aggregated MPR values after the conventional periodontal therapy in studies from Zaugg et al. (2014). These data provide a basis for applying common methods of learning from historical data (Neuenschwander et al., 2010; Neuenschwander & Schmidli, 2020; Schmidli et al., 2014).

Sample sizes reported in Table 1 are problematic, because they lead to long 95 confidence intervals (95CI) for the true mean and the true variance. For example, if the true standard deviation is fixed at 2, sample sizes lead to approximate lengths of 95CIs equal to for the true mean and for the true variance. As expected, these lengths decrease with increasing sample size. However, this decrease is faster for the mean than for the variance, indicating that the standard deviation (sd) estimates are more uncertain than the mean estimates of Table 1, especially for small total sample sizes. Consequently, the standard errors of MPR () supplied for the Bayesian meta‐analysis are also affected by this substantial uncertainty.

Binary outcome: Diabetes mellitus

Each meta‐analysis shows specific information allocation across individual studies and specific total information, which is quantified by the sum of within‐study precisions (Hardy & Thompson, 1998). Various inclusion and exclusion criteria can generate subgroups of the original data set, with some studies excluded (Bjordal et al., 2004; Crins et al., 2014; Smith et al., 1995). This can change both the total information and the information allocation across individual studies. In the Bayesian setting, the altered information allocation across individual studies can impact both sensitivity and identification estimates. However, quantification of such an impact is challenging.

To illustrate this challenge, we examine the data on the association between diabetes mellitus and the severe course of COVID‐19 (Kumar et al., 2020). To assess this association, Kumar et al. (2020) conducted a systematic review and identified observational studies, which included at least 100 adult COVID‐19 patients and reported the severity of COVID‐19. Kumar et al. (2020) determined the number of patients with the severe course of COVID‐19 in each study. They categorized the course of COVID‐19 as severe either by standard predefined criteria or when patients suffered from progressive or refractory disease, showed acute respiratory distress syndrome, required ICU care, or required invasive mechanical ventilation. Moreover, Kumar et al. (2020) determined the number of diabetic patients and computed the log‐odds ratio of severe course of COVID‐19 between diabetic and nondiabetic COVID‐19 patients and the corresponding in each study. The median number of patients included in these studies was 208 (, , , ) with a total of 13,954 COVID‐19 patients. According to Kumar et al. (2020), 23 studies had good quality, but the largest study based on registry data from 6637 patients (Study 20 in Table 2) was only classified as a fair quality study.

We follow the approach taken by Kumar et al. (2020) and consider two data sets: one based on studies (DM24, Table 2) and one based on studies (DM23) with Study 20 removed from the DM24 data set. For DM24, the total information computed as the sum of within‐study precisions of individual studies is equal to 208.1 (Hardy & Thompson, 1998). Study 20 comprises a large proportion (40.5%) of this total information. After removing Study 20, the total information of DM23 decreases to 123.9 with the maximal information contribution smaller than 14%. Considering both DM24 and DM23 data sets allows us to address the important question of how sensitivity and identification estimates are impacted by different information allocations across both data sets. However, a method is required that can quantify such an impact.

METHODS

In this section, we introduce the Bayesian NNHM and summarize the idea behind RLMC, which explicitly controls the interplay between the study‐specific within‐study standard deviation contained in the data and the impact of the heterogeneity prior in a Bayesian NNHM. In addition, we define a measure of affinity (BC) and systematic perturbations for prior () and likelihood (). Moreover, we define two combined (raw ‐ and relative ‐) sensitivity‐identification measures. Finally, we describe the implementation of the method as an R package.

Normal–normal hierarchical model

A Bayesian NNHM for studies is a Bayesian hierarchical model that consists of three parts: the sampling model (likelihood), the latent random‐effects model, and priors (Normand, 1999). The likelihood is determined by data and by the assumption that observations are realizations of a normal random variable , where and the within‐study standard deviations denoted by , for , are fixed, that is, assumed to be known. Exchangeability of the latent random effects is imposed by assuming that they follow a normal distribution with mean and a heterogeneity (between‐study) standard deviation Finally, priors for and are assumed. Given priors and , a Bayesian NNHM updates the information in the data and arrives at marginal posterior estimates for parameters .

The interplay between and in a Bayesian NNHM can be conveniently expressed by the RLMC (Ott et al., 2021), the latent model complexity (Spiegelhalter et al., 2002) per study: where the reference within‐study standard deviation is defined by the geometric mean, The RLMC for a Bayesian NNHM defined in Equation (3) attains values in independently of the primary outcome and the number of studies involved in the meta‐analysis. RLMC is a relative measure of heterogeneity that depends on both the within‐study standard deviations and a particular value of the between‐study standard deviation , which is related to (Higgins & Thompson, 2002) and shrinkage (Gelman et al., 2014; Gelman & Hill, 2007; Efron & Hastie, 2016). Moreover, RLMC is a function of both the pooling factor (Gelman & Hill, 2007) and the degree of heterogeneity (Neuenschwander et al., 2010).

The meaning of all parameters in a Bayesian NNHM in medical applications can be summarized as follows: The parameter represents the main effect. Decisions are guided by whether or not the 95 credible interval for (95CI()) contains 0. The parameter describes the between‐study heterogeneity. This is a key parameter directly impacting the model complexity of a Bayesian NNHM. The posterior distributions of latent random effects specify study‐specific effects. They have the potential to impact design and analyses of future studies when a method is applied to dynamically borrow strength from another study (Röver & Friede, 2020). The predicted true effect in a new study expresses the posterior knowledge about that study (Röver, 2020). The distribution can impact the meta‐analytic‐predictive approach to incorporating relevant historical information into the design or analysis of clinical trials (Neuenschwander & Schmidli, 2020; Schmidli et al., 2014).

To evaluate the combined evidence, we apply the Bayesian NNHM model to aggregated observations. Therefore, denotes the value of the aggregated outcome of interest in the th study and denotes its standard error. We compute Bayesian NNHMs with the R package bayesmeta (Röver, 2020).

To the prior , we assign distribution with and standard deviation , the so‐called unit‐information prior (Günhan et al., 2020; Kass & Wasserman, 1995; Röver, 2020). Note that both parameters for the prior will be kept fixed for all analyses in this paper.

For the heterogeneity prior , we assume either a light‐tailed HN or a heavy‐tailed HC distribution. The use of both HN and HC heterogeneity priors appropriately reflects the current best practice recommended for the Bayesian NNHM: HN with the base‐scale parameter fixed either at 0.5 or 1 (HN(0.5), HN(1)) (Friede et al., 2017a, 2017b; Bender et al., 2018) and HC with the base‐scale parameter fixed at 1 (HC(1)) (Gelman, 2006; Polson & Scott, 2012; Gelman et al., 2017). In the main manuscript, the choice of the base‐scale parameter for HN and HC is guided by this recommendation. In order to enable comparisons, we will align HN and HC heterogeneity priors. The base‐scale parameter values of the aligned heterogeneity priors will be obtained using an appropriate 50‐RLMC‐based adjustment of heterogeneity priors (Ott et al., 2021), reviewed in the next section. Note that the choice of base‐scale parameter values in Sections 4.1, 4.3, and 4.4 of the Supporting Information is based on the grid of 0.25, 0.5, and 0.75 RLMC values combined with the 50‐RLMC‐based adjustment.

50‐RLMC‐based adjustment of heterogeneity priors

The heterogeneity standard deviation in Equation (2) is a very important parameter for quantifying the similarity of studies. The value means that the sampling mean of all studies is (no heterogeneity and no random effects are available). In contrast, indicates an extremely large between‐study heterogeneity. In such a case, the studies are considered distinct, they must be analyzed in isolation, and different random effects are necessary to adequately model the data. Parameter is also an essential part of , the so‐called degree of heterogeneity, which is specified by the ratio of the within‐study to the between‐study standard deviations and is the key parameter determining the impact of the heterogeneity prior on the Bayesian NNHM (Neuenschwander et al., 2010).

In this paper, we will assign either an HN or an HC heterogeneity prior to and focus on the interplay between the between‐study heterogeneity prior for and the within‐study standard deviations . HN(0.5), HN(1), and HC(1) heterogeneity priors recommended in the literature induce median RLMC values 0.76, 0.93, and 0.96 for MPR (Table 3). Therefore, we anchor analyses in Section 4.1 on these three RLMC targets. In contrast, HN(0.5), HN(1), and HC(1) induce median 0.33, 0.67, and 0.82 RLMC values for DM24 (Table S16 of the Supporting Information), which we use as RLMC targets for analyses in Section 4.2. Because the HN(0.5), HN(1), and HC(1) heterogeneity priors recommended in the literature induce data‐dependent RLMC targets, we provide an alternative analysis in Sections 4.1, 4.3, and 4.4 of the Supporting Information that keeps the mutual interplay between the between‐study heterogeneity and the within‐study standard deviation values expressed by RLMC fixed on the grid of three 0.25, 0.5, and 0.75 values.

TABLE 3

Mean pocket reduction: base‐scale parameter values of 50‐RLMC‐adjusted HN and HC heterogeneity priors across three RLMC targets 0.76, 0.93, and 0.96. Targets 0.76, 0.93, and 0.96 were obtained for HN(0.5), HN(1), and HC(1) applied to MPR

Prior	RLMC=0.76	RLMC=0.93	RLMC=0.96
HN	0.50	1.00	1.48
HC	0.34	0.67	1.00

Given a target RLMC value, Ott et al. (2021) introduced a 50‐RLMC‐based adjustment for heterogeneity priors, which aligns the medians of different heterogeneity priors at the reference threshold . Thus, where the reference threshold satisfies Note that the condition for the tail adjustment in Equation (5) aligns the medians of different heterogeneity priors at a dynamic reference threshold in Equation (6). is dynamic, because it depends on both the data through and the base target RLMC. Because the densities of HN and HC for are two different mathematical functions, it is impossible to get a perfect match. Nevertheless, the 50‐RLMC‐based adjustment uses to align the medians of both the HN and HC heterogeneity priors.

Given the above prerequisites, the heterogeneity prior base‐scale parameter values in Table 3 and in Table S16 of the Supporting Information are found numerically to fit the condition imposed by Equation (5). For any data set at hand, there is a one‐to‐one correspondence between the base‐scale parameter value of the heterogeneity prior and the target RLMC value. For HN and HC, formulas connecting the tail constraint in Equation (5) and the base‐scale parameter value of the heterogeneity prior can be derived analytically to get for HN and for HC (Ott et al., 2021). For example, if we assume base for the MPR data with in Table 1, the scale parameter for a 50‐RLMC‐adjusted HN is 0.5 (Table 3). For the base and the same MPR data set, the scale parameter for a 50‐RLMC‐adjusted HC is specified as 0.34. This approach enables a coherent choice of the scale parameter values of two different HN and HC heterogeneity priors (Ott et al., 2021). Our experience so far reveals that it is more meaningful to compare posteriors induced by 50%‐RLMC‐adjusted heterogeneity priors, which represent a similar model complexity. Therefore, Figure 1 is split into panels that show posteriors across different RLMC targets.

FIGURE 1

Mean pocket reduction: prior densities for (top row) and posterior densities for (middle row) and (bottom row) induced by 50‐RLMC‐adjusted HN (solid line) and HC (dashed line) heterogeneity priors with target RLMC values fixed at 0.76 (left column), 0.93 (middle column), and 0.96 (right column). See Table 3 for the corresponding base‐scale parameter values of HN and HC heterogeneity priors

Affinity measure BC

To quantify the impact of data and heterogeneity prior perturbations on the posterior inference in a Bayesian NNHM, we need a formal measure to assess the affinity of two distributions. We suggest using the Bhattacharyya coefficient (BC), which is a symmetric measure of affinity (Roos et al., 2015; Roos & Held, 2011). The BC between a base (b) and an altered (a) density function is given by with . BC is equal to the maximal value 1 if and only if the densities are identical. It takes the minimal value 0 when the densities do not overlap due to disjoint supports.

BC has convenient numerical properties and is invariant to any one‐to‐one transformation of both densities (Roos & Held, 2011; Jeffreys, 1961). Particularly significant is that the value of BC computed for two densities and two log‐transformed densities is equal. BC is also the main component of the Hellinger distance (H), which naturally fits into the framework of Bayesian geometry (de Carvalho et al., 2019). More specifically, .

The BC in Equation (7) between two normal densities can be derived analytically. Denote by the density of a normal distribution, for . Then, the total BC for two normal densities is

When computing BC in practice, we utilize the general idea that any density can be approximated to the first order by a parametric normal distribution (Johnson et al., 1994; Hjort & Glad, 1995). We thus focus our attention on the parametric normal approximation of the marginal posterior distribution and approximate for by an appropriate normal distribution combined with moment matching and use Equation (8) to get Moment matching equates the mean and standard deviation estimates of a marginal posterior distribution and the approximating normal distribution . Note that for the parameter , we apply the normal approximation to the log‐transformed marginal posterior for . This approach is justified by the invariance of BC to any one‐to‐one transformation. For brevity of notation, we refer to this parameter as , although all numerical computations are conducted based on .

Note that Equation (8) applied to in Equation (9) focuses on changes in posterior moments of base and altered marginal posterior distributions. In fact, the total BC in Equation (8) is based on common descriptive statistics summarizing posterior inference: posterior mean and posterior standard deviation. Equation (8) quantifies the impact of perturbations on both posterior moments in two steps. While the multiplicator measures the discrepancy between two standard deviations, the multiplicand measures the discrepancy between two means that is scaled by the pooled variance. Note that the multiplicand is the well‐known Mahalanobis distance in one dimension (Mardia et al., 1979).

Systematic heterogeneity prior and likelihood perturbations

Sensitivity () assesses to what extent systematic heterogeneity prior perturbations cause posterior variation. In contrast, identification () reflects how perturbations of the likelihood () impact posterior inference, which can be achieved by suitably scaling the values in .

The posterior distribution underlying the Bayesian meta‐analysis expressed by an NNHM is a self‐contained system in Equation (10). This system consists of the likelihood (), which encompasses the data, the latent random‐effects model and priors () (see Section 2.1 of the Supporting Information for details): Note that the within‐study standard deviation values in the data specify the impact of the likelihood defined in Equation (1) on the posterior in Equation (10). The heterogeneity prior specification governs the impact of the heterogeneity prior on the posterior in Equation (10).

The main objective of our method for a Bayesian NNHM is to coherently assess the impact of two different sources of uncertainty: the uncertainty of heterogeneity prior choice and the uncertainty in within‐study standard deviation values. The impact of both sources is assessed by appropriately chosen systematic perturbations. In the first dimension (), the uncertainty inherent in a specific choice of the scale parameter value of a heterogeneity prior is assessed using appropriate scale parameter value perturbations. In the second dimension (), the uncertainty inherent in within‐study standard deviation values provided by the data is assessed using appropriate value perturbations. Moreover, perturbations in both dimensions are unified and made coherent, because the impact of a heterogeneity prior for and the within‐study standard deviation values is connected through the RLMC in a Bayesian NNHM (see Section 2.2 of the Supporting Information for details). For example, lower values of RLMC in Equation (3) can be obtained either by perturbing the heterogeneity prior () with fixed data or perturbing data in the likelihood () with a fixed heterogeneity prior. This section describes systematic perturbations of heterogeneity priors () and the likelihood ().

Systematic heterogeneity prior () perturbation

Systematic perturbations of the heterogeneity prior () in a Bayesian NNHM are conducted given a fixed likelihood. In particular, both the within‐study standard deviation values in the data and consequently the reference within‐study standard deviation in Equation (4) are kept fixed. If RLMC changes while the likelihood is kept fixed, then the scale parameter of a heterogeneity prior for must change. Systematic perturbations of the heterogeneity prior () are carried out using the 50‐RLMC‐based adjustment. Perturbing, that is, slightly changing the value of the RLMC, leads to slight changes of the reference threshold according to Equation (6). Because the tail adjustment with the tail fixed at in Equation (5) is directly connected through an analytical formula to the scale parameter of heterogeneity priors, slight perturbations of imply slight perturbations of the scale parameter of heterogeneity priors.

As an example, we demonstrate the method described above on the MPR data in Table 1 and show how the scale parameter of an HN heterogeneity prior is perturbed through an RLMC perturbation with fixed values. If we assume the base , then we obtain the base‐scale parameter value 0.502 for an HN heterogeneity prior: HN(0.502). In addition, with , we obtain , which implies HN(0.496), and , which implies HN(0.508).

Systematic likelihood () perturbation

Here, we show how systematic perturbations of RLMC provoke systematic perturbations of the likelihood () (data), given a fixed RLMC‐adjusted heterogeneity prior. In order to perturb the likelihood in a Bayesian NNHM subject to the RLMC constraint, we keep the heterogeneity prior fixed according to the base 50‐RLMC‐adjustment and multiply the within‐study standard deviations by a factor . Thus, where with altered values: either or . Scaling the data by a factor according to Equations (11) and (12) while keeping the heterogeneity prior fixed at a base target makes the likelihood more or less pronounced. This corresponds to changing the effective sample size in each study supplied for a Bayesian meta‐analysis.

For example, if we set the base and take , then we obtain the lower , which implies a factor and the upper , which implies a factor . Changing RLMC with a fixed heterogeneity prior will induce scaling of the study‐specific within‐study standard deviation values by the factor , and data with will be supplied for a Bayesian meta‐analysis.

Both perturbed heterogeneity priors and perturbed study‐specific within‐study standard deviation values used for Bayesian NNHM computations produce perturbed marginal posterior distributions of all parameters in the Bayesian NNHM . Therefore, we are able to compute BC between base and altered marginal posteriors of with respect to the RLMC perturbations affecting either the heterogeneity prior or the within‐study standard deviation . This property leads to the unified two‐dimensional sensitivity‐identification (‐) measure, which is defined in the next section as a numerical approximation to the second derivative of BC with respect to systematic RLMC perturbations evaluated at a base target RLMC. Note that RLMC targets are fixed at 0.76, 0.93, and 0.96 for MPR in Table 3, at 0.33, 0.67, and 0.82 for DM in Table S16 of the Supporting Information, and at 0.25, 0.5, and 0.75 for analyses reported in Sections 4.1, 4.3, and 4.4 of the Supporting Information.

Formal combined sensitivity‐identification (‐) measure

Formal approaches to sensitivity and identification quantification are usually based on differential calculus or its numerical approximations (McCulloch, 1989; Xie & Carlin, 2006). For the BC, the second derivative is a practical quantity (McCulloch, 1989; Dey & Birmiwal, 1994).

We define the unified, formal, two‐dimensional sensitivity‐identification (‐) measure as negative second derivatives of the BC with respect to the RLMC perturbations evaluated at a base value when perturbing either the heterogeneity prior () or the , , values in the likelihood () for . The negative second derivatives of BC in Equations (13) and (14) with respect to RLMC perturbations assess the curvature of the squared Hellinger distance (). They quantify the acceleration with which the marginal posterior changes locally around when RLMC perturbations are induced either by or by perturbations. We compute numerical approximations to the derivatives and approximate the posterior distribution with a normal distribution when computing BC as shown in Equation (9). Thus, we substitute the true second derivative with its numerical approximation Note that the and estimates are always positive, because BC attains a maximum 1 at . Large values of both ‐ estimates indicate that marginal posteriors are highly affected both by the uncertainty in the choice of the heterogeneity prior scale parameter ( dimension) and by the uncertainty of the within‐study standard deviations ( dimension).

To compute the raw and estimates in practice, we compute the Bayesian NNHM several times. First, for a fixed base value, the base Bayesian NNHM for the original base data with a 50‐RLMC‐adjusted heterogeneity prior is computed. This leads to base marginal posterior distributions for all parameters . Second, we systematically perturb the RLMC and compute the posteriors that are induced, either when the systematic RLMC perturbation impacts the likelihood (through perturbed , values) with a fixed heterogeneity prior or when it impacts the heterogeneity prior with fixed data. Note that the value of the affinity measure for two nonperturbed posteriors.

We then apply the combined ‐ measure, anchored at a base RLMC, to both HN and HC heterogeneity priors. For numerical computations in this paper, we perturb each RLMC target by a value . For example, if we set the base target and perturb it by the perturbation 0.0044, then we obtain two perturbed RLMC values: and . We postpone the justification of this particular value choice and a discussion of the stability of the resulting ‐ values until Section 6 of the Supporting Information. Although the analytical second derivative has been shown to be a well‐defined and practical quantity (McCulloch, 1989; Dey & Birmiwal, 1994), its numerical approximation can potentially attain a value of 0. However, this issue did not turn out to be a problem in the medical applications considered in this paper.

Quotients of sensitivity identification (‐)

The unified ‐ measure defined in Section 3.7 allows for a wide range of direct comparisons relevant for applications. Such comparisons are expressed by quotients of ‐ values leading to relative ‐ estimates. For example, one set of relevant ‐ values emerges when we wish to compare the impact of two different heterogeneity priors on each parameter in the Bayesian NNHM. In such a case, the unified two‐dimensional relative ‐ measure is defined by with and a fixed RLMC target.

Another interesting ‐ measure is obtained when raw ‐ values are compared between two different data sets. For example, ‐ values for obtained for a data set based on studies can be compared with those obtained for a subgroup based on studies. Thus, the formal unified two‐dimensional relative ‐ measure is defined by with for a fixed target RLMC and a fixed heterogeneity prior.

Depending on the research question, other relative ‐ estimates can be computed. For example, it can be relevant to compare the impact of two different Bayesian NNHM and NtHM models for meta‐analysis or the impact of randomized and nonrandomized subgroups (Section 2.6 of the Supporting Information). Alternatively, it can be of interest to compare the impact of two different RLMC targets. Moreover, raw ‐ values for each parameter can be compared within one scenario, with both the target RLMC and the heterogeneity prior kept fixed, as defined in Section 2.6 of the Supporting Information.

Package: si4bayesmeta

The method described above is freely accessible as an R package si4bayesmeta. The si4bayesmeta package (https://github.com/hunansona/si4bayesmeta) contains several wrapper functions built on two R packages bayesmeta (Röver, 2020) and rjags (Plummer, 2016). Functions median, pri, and pri_par compute RLMC targets and base‐scale parameters of HN and HC heterogeneity priors shown in Table 3 and in Table S16 of the Supporting Information. Two functions, d2BC_S_I_HN_raw and d2BC, produce the main result for Bayesian NNHM fitted by bayesmeta shown in Table 5 below and in Table S18 of Supporting Information for a specified RLMC target. Posterior estimates for a specified RLMC target in Table 4 below and in Table S17 of the Supporting Information are computed by the functions raw and raw. Additional functions necessary to compute ‐ estimates for Bayesian NNHM and Bayesian NtHM fitted by MCMC sampling in JAGS are specified in Section 2.5 of the Supporting Information. Moreover, the function h facilitates the choice of the numerical RLMC perturbation , which is adjusted for the epsilon‐local sensitivity (Roos et al., 2015), as outlined in Section 6 of the Supporting Information. All functions are applicable to data frames in the bayesmeta format. Source code to reproduce the results is available as Supporting Information on the journal's web page. (http://onlinelibrary.wiley.com/doi/xxx/suppinfo)

TABLE 5

Mean pocket reduction ‐: sensitivity and identification estimates for parameters in Bayesian NNHMs obtained for six scenarios across target RLMC values fixed at 0.76, 0.93, and 0.96 with HN (upper batch) and HC (lower batch) 50‐RLMC‐adjusted heterogeneity priors. See Table 3 for the corresponding base‐scale parameter values of HN and HC heterogeneity priors

	RLMC=0.76		RLMC=0.93		RLMC=0.96
Par	S_P	I_L	S_P	I_L	S_P	I_L
μ^(HN)	0.1661	1.3805	0.2102	6.7392	0.2486	24.8017
τ^(HN)	3.4663	470.1564	2.3836	3393.3009	70.0537	7766.1207
θnew(HN)	0.1855	0.2790	0.2458	1.5349	0.2940	5.8388
μ^(HC)	0.0039	1.2579	0.0800	7.6894	0.3580	28.1756
τ^(HC)	0.5671	307.8816	0.1989	1099.9935	5.0586	2276.8277
θnew(HC)	0.0043	0.3043	0.0903	1.8186	0.4139	6.5216

TABLE 4

Mean pocket reduction: posterior mean and standard deviation (sd) estimates together with the shortest 95 credible interval (CI) and the length it reaches (LCI) for parameters obtained with target RLMC values fixed at 0.76, 0.93, and 0.96 and with an HN (upper batch) and an HC (lower batch) 50‐RLMC‐adjusted heterogeneity prior. See Table 3 for the corresponding base‐scale parameter values of HN and HC heterogeneity priors

	RLMC=0.76					RLMC=0.93					RLMC=0.96
Par	Mean	sd	CI_low	CI_up	LCI	Mean	sd	CI_low	CI_up	LCI	Mean	sd	CI_low	CI_up	LCI
μ^(HN)	1.14	0.22	0.71	1.58	0.87	1.15	0.24	0.68	1.64	0.96	1.15	0.25	0.67	1.65	0.99
τ^(HN)	0.71	0.16	0.43	1.02	0.60	0.79	0.20	0.44	1.18	0.74	0.81	0.21	0.44	1.23	0.79
θnew(HN)	1.14	0.76	−0.37	2.66	3.03	1.15	0.85	−0.53	2.85	3.38	1.15	0.87	−0.58	2.91	3.48
μ^(HC)	1.14	0.23	0.69	1.61	0.91	1.15	0.24	0.69	1.62	0.94	1.15	0.24	0.68	1.64	0.96
τ^(HC)	0.74	0.19	0.42	1.12	0.70	0.76	0.20	0.43	1.16	0.73	0.78	0.20	0.44	1.19	0.75
θnew(HC)	1.14	0.80	−0.44	2.75	3.19	1.15	0.82	−0.49	2.81	3.30	1.15	0.84	−0.53	2.85	3.38

RESULTS

In this section, we provide posterior inference, ‐ values, and relative ‐ estimates for parameters for MPR (Table 1) and for both DM24 and DM23 in Table 2. Here, we consider six scenarios that arise for three RLMC targets induced by HN(0.5), HN(1), and HC(1) heterogeneity priors crossed with two 50‐RLMC‐tail adjusted heterogeneity priors (HN, HC). The Supporting Information extends the results of the main manuscript and considers several additional topics, which are reviewed in Section 4.3.

Mean pocket reduction: Comparison of HN and HC heterogeneity priors

For historical MPR data (, Table 1), HN(0.5), HN(1), and HC(1) heterogeneity priors, which are recommended in the literature, induce the RLMC equal to 0.76, 0.93, and 0.96, respectively. Table 3 focuses on these target RLMC values and reports the base‐scale parameter values of the 50‐RLMC‐adjusted HN and HC heterogeneity priors. The densities of these HN and HC priors for are shown in the top row of Figure 1. It is evident that the spread of both HN and HC increases with increasing base‐scale parameter values and target RLMC values. The resulting posterior densities of (middle row of Figure 1) and (bottom row of Figure 1) induced by the HN and HC heterogeneity priors agree well across target 0.76, 0.93, and 0.96 RLMC values.

Table 4 complements these visual findings and provides posterior mean and standard deviation estimates for parameters together with the shortest 95 credible interval and the length we obtained for it in the six different scenarios. For both HN and HC, the length of 95CIs for the main effect , the heterogeneity , and the predicted true effect in a new study increase with increasing RLMC. The 95CIs are longer for HC than for HN for but are shorter across the 0.93 and 0.96 RLMC targets.

For the main effect and both the 0.76 and 0.93 RLMC targets, HC leads to low values in Table 5. This property is reversed when . The marginal posterior for is more affected by HN than by HC heterogeneity prior perturbations.

Note that for all 0.76, 0.93, and 0.96 RLMC targets, the marginal posterior for is highly affected by the ‐uncertainty for both HN and HC heterogeneity priors (Table 5). Moreover, the posteriors of all parameters are increasingly affected by ‐uncertainty as RLMC values increase.

To assess how the impact of the HN and HC heterogeneity priors differs, Equation (15) has been applied to ‐ values of Table 5, leading to relative ‐ estimates in Table 6. For , ‐ estimates are greater than 1, revealing that the impact of both the heterogeneity prior and the within‐study standard deviation uncertainty on parameters is less pronounced for HC than for the HN heterogeneity prior. This relation for and is reversed for with ‐ estimates less than 1. For , ‐ estimates show that posterior estimates of and obtained from HN are 2.63 and 2.72 times more sensitive than HC but they are less affected (0.88, 0.84) by standard deviation uncertainty. This indicates that for HN(1) and HC(0.67) corresponding to , the advantages of the HN and the HC heterogeneity priors balance out and the decision as to which of the heterogeneity priors should be recommended for analysis of MPR is less clear‐cut.

TABLE 6

Mean pocket reduction ‐: quotients of sensitivity and identification estimates for HN with respect to HC defined in Equation (15) for parameters in Bayesian NNHMs across three RLMC targets 0.76, 0.93, and 0.96 with 50‐RLMC‐adjusted HN and HC heterogeneity priors. See Table 3 for the corresponding base‐scale parameter values of HN and HC heterogeneity priors

	RLMC=0.76		RLMC=0.93		RLMC=0.96
Par	qS^HN/HC	qI^HN/HC	qS^HN/HC	qI^HN/HC	qS^HN/HC	qI^HN/HC
μ	42.51	1.10	2.63	0.88	0.69	0.88
τ	6.11	1.53	11.99	3.08	13.85	3.41
θ_new	43.08	0.92	2.72	0.84	0.71	0.90

The ‐ approach can be used to explore other RLMC targets such as, for example, 0.25, 0.5, and 0.75 in Section 4.1 of the Supporting Information. Note that Table 5 of the main manuscript and Table S9 of the Supporting Information together cover the range of RLMC targets between 0.25 and 0.96 and the range of base‐scale parameters from 0.16 to 1.48 for HN and from 0.11 to 1 for HC. They show that HC is preferable for MPR data for a wide range of RLMC targets.

One advantage of the unified ‐ approach is that it clearly demonstrates that depending on the RLMC value, the posterior estimates are affected to different extent by the heterogeneity prior uncertainty and the within‐study standard deviation uncertainty. Moreover, it provides some guidance as to which Bayesian model specifications are beneficial. This can be particularly useful in cases where substantial information about the between‐study heterogeneity () is missing. For example, a practitioner concerned with the small sample sizes of individual studies and the validity of the corresponding within‐study standard deviation estimates of historical MPR data (Table 1) could choose an HC(0.19) () heterogeneity prior for Bayesian analysis, because it renders the lowest estimate for in Table S9 of the Supporting Information. On the other hand, when a practitioner prefers to use the HC(1) heterogeneity prior recommended in the literature ( in Table 3), he/she should be aware that the validity of the within‐study standard deviations in Table 1 is very important for this choice, because the posterior of is highly affected by the uncertainty in within‐study standard deviation assumptions (Table 5). Moreover, he/she can now be aware that the use of an HN(1.48) rather than HC(1) heterogeneity prior would be recommended for this target (Table 5).

Diabetes mellitus: Comparison of 24 and 23 studies

HN(0.5), HN(1), and HC(1) applied to DM24 data in Table 2 induced median RLMC values equal to 0.33, 0.67, and 0.82 (Table S16 of the Supporting Information). Note that these RLMC values differ from the RLMC values obtained for the same HN(0.5), HN(1), and HC(1) heterogeneity priors for MPR data analyzed in Section 4.1. For DM24 and DM23, we fixed target RLMC values at 0.33, 0.67, and 0.82 and adjusted the base‐scale parameters of HN and HC heterogeneity priors with the 50‐RLMC‐based adjustment (Table S16 of the Supporting Information). There is a one‐to‐one correspondence between base scale parameters of HN and HC heterogeneity priors and RLMC targets. Larger values are linked to larger RLMC and lead to prior density distributions for that demonstrate larger spread in Figures 2 and 3 (top row) of the Supporting Information. The resulting posterior distributions for (middle row) and (bottom row) agree well across 0.33, 0.67, and 0.82 RLMC targets. This visual impression is confirmed by stable posterior estimates reported in Table S17 of the Supporting Information.

Table S18 of the Supporting Information shows ‐ estimates for both DM24 (upper batch) and DM23 (lower batch) data sets. Although both sample sizes and 24 are large, the posterior for is sensitive to the within‐study standard deviation uncertainty for both HN and HC heterogeneity priors. Removal of the large Study 20 leads to larger ‐ estimates for DM23 than those for DM24, as indicated by ‐ estimates in Table S19 of the Supporting Information. This removal affects both and dimensions for the overall mean parameter but the dimension is affected more. Posteriors of all parameters for DM23 are more affected by within‐study standard deviation uncertainty than the parameters for DM24. In particular, the removal of Study 20 leads to posteriors for that are more than eight times more affected by within‐study standard deviation uncertainty than posteriors for DM24. The sensitivity to heterogeneity prior uncertainty of posteriors induced by HC is more affected by the removal of Study 20 than the sensitivity of posteriors induced by HN.

For DM24, HC heterogeneity prior is of advantage only for the RLMC target equal to 0.33 (see ‐ estimates in Table S20 of the Supporting Information). Otherwise, the use of the HN heterogeneity prior induces posteriors that are less affected by the heterogeneity prior uncertainty and the within‐study standard deviation uncertainty for both DM24 and DM23. Note that Tables S23– S25 of the Supporting Information based on 0.25, 0.5, and 0.75 RLMC targets provide similar results.

These results complement the analysis by Kumar et al. (2020) for both DM24 and DM23 data sets and equip them with ‐ estimates. One advantage of the unified ‐ approach is that we can explicitly quantify how the change of information allocation across individual studies after removal of one large study impacts the sensitivity of posteriors to both the heterogeneity prior uncertainty and the within‐study standard deviation uncertainty on an equal footing. Moreover, if no substantial information about the between‐study heterogeneity () is available, we can recommend the use of HN(1) () and HN(1.48) () for the analysis of both DM24 and DM23 data sets.

Summary of results provided in the Supporting Information

The Supporting Information extends the results of the main manuscript in several ways. First, it focuses on the grid of 0.25, 0.5, and 0.75 RLMC targets and demonstrates that the ‐ approach can be anchored either on heterogeneity prior values known from the literature (main manuscript) or directly on RLMC targets. Numerical issues related to the choice of the RLMC perturbation and the stability of estimates are addressed in Section 6 of the Supporting Information.

Second, Section 2.7 of the Supporting Information compares ‐ with two other approaches to sensitivity assessment: the informal approach (Gustafson, 1996) and the epsilon‐local sensitivity measure (Roos et al., 2015). While the informal approach has no ability to quantify and compare sensitivity between different specifications of Bayesian meta‐analyses, the epsilon‐local sensitivity measure focuses on data‐to‐no‐data comparison and is therefore restricted to two parameters ( and ). In contrast, ‐ quantifies the impact of base‐scale parameter uncertainty of HN and HC heterogeneity priors and within‐study standard deviation uncertainty on posterior estimates of all parameters in Bayesian meta‐analysis.

Third, Section 4.4 of the Supporting Information provides ‐ estimates for one subgroup AGR2 based on randomized studies and another subgroup AGR4 based on nonrandomized studies. These ‐ estimates warn that properties of randomized AGR2 and nonrandomized AGR4 data sets differ. This indicates that a more sophisticated model than the simple Bayesian NNHM would be better suited for the combined analysis of both disparate types of evidence. For example, the Bayesian generalized evidence synthesis (GES) based on a three‐level hierarchical model, which is designed to synthesize disparate types of evidence, would be preferable for the analysis of the combined data set based on six observations (Prevost et al., 2000; Sutton & Abrams, 2001; Verde & Thomas, 2015; Verde, 2021).

Fourth, Sections 2.4 and 2.5 of the Supporting Information extend the ‐ approach to Bayesian NtHM, which is useful when the assumption of the normality of random effects is violated. Section 4.1 of the Supporting Information extends the analysis of MPR data from Section 4.1 of the main manuscript and investigates whether using Bayesian NtHM has advantages over using Bayesian NNHM for this data set. Table S12 of the Supporting Information demonstrates that applying Bayesian NtHM to MPR data produces posterior estimates that are less affected by both the heterogeneity prior uncertainty and the within‐study uncertainty than posteriors produced by Bayesian NNHM across 0.25, 0.5, and 0.75 RLMC targets. Table S14 (lower part) confirms that for 0.25, 0.5, and 0.75 RLMC targets, the advantage of using HC rather than HN also applies to the Bayesian NtHM.

Finally, Section 5 of the Supporting Information extends one simulation design suggested by Hardy and Thompson (1998) to cover the range of cases, when the assumption of the normality of random effects is violated. This simulation study based on 1000 simulations for each scenario shows that ‐ can correctly demonstrate that the Bayesian NtHM leads to lower ‐ estimates than the Bayesian NNHM for both HN and HC heterogeneity priors.

DISCUSSION

We showed that HC was preferable when summarizing the historical MPR data for the conventional periodontal therapy. These results demonstrate that historical data carry more information than what can be summarized by plain estimates, knowledge that provides some guidance on the choice of the heterogeneity prior for further analysis. The extensive use of all of the information provided by historical data is particularly relevant for areas of research for which obtaining patients is costly and where meta‐analyses commonly contain small sample sizes. Thus, ‐ estimates could support the application of common methods of learning from historical data (Neuenschwander et al., 2010; Schmidli et al., 2014; Neuenschwander & Schmidli, 2020), in evidence‐based dentistry and beyond. More work is, however, necessary to determine if ‐ estimates could also help judge whether mean difference or standardized mean difference is better suited to analyze continuous outcomes in applications (Hedges, 1981; Sedgwick & Marston, 2013; White & Thomas, 2005).

The proposed two‐dimensional ‐ system for sensitivity comparisons is designed for the conventional NNHM model for aggregated data, which is widely used in applications (Normand, 1999; Gelman et al., 2014; Hoaglin, 2015; Jackson et al., 2018). Note that the accuracy of ‐ estimates depends on how well the NNHM assumptions match the data at hand, and this match is expected to be better for continuous outcomes than for noncontinuous ones. Although Spiegelhalter et al. (2004) argue that the normal approximation to log‐OR scale is reasonable, it should be noted that ‐ values for binary DM and AGR data are, in fact, an approximation.

We considered the DM24 data set dealing with the association of diabetes mellitus and the severe course of COVID‐19, explicitly quantifying the impact on ‐ estimates of omitting one large study in DM23. Moreover, we showed that applying an HN heterogeneity prior for the analysis of DM24 and DM23 has the advantage of producing posterior estimates that are less affected both by the heterogeneity prior and by the within‐study uncertainty. Although the advantage of using HN rather than HC was clear for DM24 and DM23, HC was preferable when analyzing MPR. These results show that the HC heterogeneity prior (Gelman, 2006; Polson & Scott, 2012; Gelman et al., 2017) is a valid alternative to HN, which is often recommended in the literature (Bender et al., 2018; Friede et al., 2017a, 2017b). Future work aims to determine under which specific conditions HC and HN are beneficial in applications.

At first glance, the method for taking derivatives of the BC with respect to RLMC perturbations, resulting in a unified, two‐dimensional ‐ measure for the Bayesian NNHM, may seem complicated. However, we are aiming for a broadly applicable, unified, two‐dimensional method, independent of the scale of the primary outcome, and achieving such general applicability necessitates the use of these mathematical tools. Taking derivatives of a distance between two posterior distributions is a well‐known tool used in Bayesian robustness (Dey & Birmiwal, 1994; Gustafson, 2000; McCulloch, 1989). The particular choice of the BC distance is motivated by its boundedness, which renders better numerical properties than those of other known distances, such as the unbounded Kullback–Leibler distance (Jeffreys, 1961; Roos & Held, 2011; Roos et al., 2015). RLMC transforms the relation between the between‐study heterogeneity and the within‐study standard deviations into a dimensionless bounded value, independent of the scale of the primary outcome. Moreover, RLMC controls the flow of mass between the between‐study heterogeneity and the within‐study standard deviations. These properties make RLMC a valid basis for the unified, two‐dimensional ‐ measure that assesses the impact of both the heterogeneity prior and the within‐study standard deviations on the posterior results.

The unified ‐ approach proposed is similar in spirit to the approach taken by Hardy and Thompson (1998) and Gustafson et al. (2006), who show that the between‐study heterogeneity can only be well understood when the within‐study standard deviation values are taken into account. In our opinion, the mutual interplay between the heterogeneity and the within‐study standard deviation cannot be easily disentangled, so we choose to simultaneously quantify the impact of and . However, the nature of the sensitivity analysis (based on prior perturbations) is very different from that of the identification analysis (based on perturbation of within‐study standard deviation values) . While the sensitivity analysis can guide in specifying the more appropriate prior on the heterogeneity parameter, the within‐study standard deviation values are provided by the data and cannot be modified. Nonetheless, the identification analysis is important, because it serves as a reminder that the within‐study standard deviation values are not true values but rather uncertain estimates, and this uncertainty can also affect posteriors.

One advantage of the unified ‐ approach is that ‐ estimates facilitate direct comparisons of sensitivity estimates induced by different model specifications on an equal footing. ‐ estimates can be easily transformed into relative ‐ estimates based on quotients. These quotients compare two specifications of the Bayesian meta‐analysis and determine which of them leads to posteriors that are less affected by uncertainty in the prior choice and the uncertainty of the within‐study standard deviation values. These comparisons may include different heterogeneity prior distributions (HN and HC), data with different allocation of information across individual studies after removal of one large study ( and ), different latent model assumptions (NNHM and NtHM), and subgroups of randomized and nonrandomized studies. Other comparisons can be designed depending on the scientific question of interest.

Although NNHM is popular among researchers, the assumption of the normality of random effects is sometimes violated in practice. Therefore, we extended the unified ‐ approach to cover Bayesian NtHM fitted by MCMC in JAGS. To compute the corresponding ‐ estimates, JAGS runs models several times based on the number of chains, the length of burn‐in, and thinning prespecified by the user and does not conduct any convergence analysis. Because these computations can be time‐consuming, more work is necessary to speed up procedures based on MCMC samples. Another approach to reduce the time necessary for computation of ‐ estimates in Bayesian NtHM would be the implementation of the Bayesian NtHM in INLA (Martins & Rue, 2014; Martins et al., 2013; Rue et al., 2009, 2017). This approach is promising because procedures based on Bayesian numerical approximation, such as, for example, bayesmeta for Bayesian NNHM, require less computational time and postprocessing, and fewer specifications to fit the model.

The unified ‐ approach could be further extended in at least two ways. One promising extension would be the application of the ‐ approach to the Bayesian GES based on a three‐level hierarchical model fitted by MCMC (Prevost et al., 2000; Sutton & Abrams, 2001). This Bayesian GES considers one between‐study heterogeneity, two type‐specific between‐study heterogeneities, and three heterogeneity priors. Although a generalized RLMC can be derived and systematic perturbations of data (‐dimension) and three heterogeneity priors (‐dimension) can be designed for this Bayesian GES, more work is necessary to implement these systematic perturbations in R and to provide the full ‐ functionality for users of this Bayesian GES. A successful implementation would provide a basis for further extensions to the modern Bayesian bias‐corrected meta‐analysis (Verde, 2021).

Another promising extension would be the use of study‐specific factors for within‐study standard deviation perturbation. Currently, likelihood perturbations are tuned by a factor , which produces the same percent change on each within‐study standard deviation. In applications, different within‐study standard deviation values are quite common, such as, for example, when evidence of many small studies is combined with a few large studies (Hardy & Thompson, 1998). To better reflect such cases, our method could be extended to allow more complex perturbations that tune various within‐study standard deviation values differently. However, allocation of different tuning factors to each study would be very complicated, because endless different combinations of tuning could be considered (Hardy & Thompson, 1998). In particular, general rules should be developed for choosing studies that require such specific tuning. This could be done based on a threshold of, for example, 50% contribution of a specific study to the total information (Hardy & Thompson, 1998) or based on the type of evidence, such as randomized or nonrandomized study (Prevost et al., 2000; Sutton & Abrams, 2001). More work is necessary to implement this new idea in its full generality in R.

The method proposed is provided as an R package si4bayesmeta openly accessible on R‐Forge. We believe that this software will facilitate the application of the unified approach to sensitivity and identification quantification in a Bayesian meta‐analysis in a wide range of contexts and in many research fields.

CONFLICT OF INTEREST

The authors have declared no conflicts of interest.

OPEN RESEARCH BADGES

This article has earned an Open Data badge for making publicly available the digitally‐shareable data necessary to reproduce the reported results. The data is available in the Supporting Information section.

This article has earned an open data badge “Reproducible Research” for making publicly available the code necessary to reproduce the reported results. The results reported in this article could fully be reproduced.

Supporting Information for this article is available on the WWW under https://doi.org/10.1002/bimj.202000193

Click here for additional data file.

Supporting information

Click here for additional data file.