Home Bayesian Variable Selection Methods for Matched Case-Control Studies
Article Publicly Available

Bayesian Variable Selection Methods for Matched Case-Control Studies

  • Josephine Asafu-Adjei EMAIL logo , Mahlet G. Tadesse , Brent Coull , Raji Balasubramanian , Michael Lev , Lee Schwamm and Rebecca Betensky
Published/Copyright: January 31, 2017
Download Article (PDF)
Become an author with De Gruyter Brill
Submit Manuscript Author Information
The International Journal of Biostatistics
From the journal The International Journal of Biostatistics Volume 13 Issue 1

Abstract

Matched case-control designs are currently used in many biomedical applications. To ensure high efficiency and statistical power in identifying features that best discriminate cases from controls, it is important to account for the use of matched designs. However, in the setting of high dimensional data, few variable selection methods account for matching. Bayesian approaches to variable selection have several advantages, including the fact that such approaches visit a wider range of model subsets. In this paper, we propose a variable selection method to account for case-control matching in a Bayesian context and apply it using simulation studies, a matched brain imaging study conducted at Massachusetts General Hospital, and a matched cardiovascular biomarker study conducted by the High Risk Plaque Initiative.

Keywords: Bayesian analysis; conditional logistic regression; matched case-control studies; variable selection methods

1 Introduction

In matched case-control studies, subjects from a particular diagnostic group(s) (cases) are matched with those from a comparison group(s) (controls) based on important demographic characteristics, such as age or gender. Matching on these potential confounders can result in substantial improvements in efficiency and statistical power to identify feature variables that are associated with case-control status [1]. Matched case-control studies are becoming increasingly popular in studies involving high dimensional data that aim to identify subsets of relevant features. Examples include brain imaging studies aimed at identifying brain regions associated with comorbidities or genomic studies focused on discovery of cancer biomarkers. Despite the popularity of matched high dimensional studies, it is quite common for these studies to ignore the matched design used when applying variable selection techniques (e.g., Anglim et al. [2], Westman et al. [3]). Failure to account for matching has been shown to decrease variable selection accuracy [1] and lead to biased results [4].

Currently, there are several frequentist variable selection approaches for matched high dimensional data that incorporate case-control matching. Tan et al. [5] develop a modified paired t-test statistic to identify a subset of relevant features that serves as a basis for classification via support vector machines (SVM). Although matching is accounted for with respect to variable selection, it is ignored with respect to building the SVM classifier. In addition, their approach of identifying relevant features involves univariate tests, which do not control for the effects of other features and, thus, can lead to spurious identification of relevant features. Adewale et al. [6] develop two modified versions of boosting for correlated binary response data. The first version utilizes a loss function for the generic gradient descent boosting algorithm [7] that handles correlated binary responses. The second version modifies the likelihood optimization boosting algorithm [7] via a generalized linear mixed modeling approach in order to handle correlated binary responses. However, boosting approaches may have trouble identifying interactions among different features [1] and have decreased accuracy for data sets with relatively small sample sizes [8]. In matched studies, a standard analytic approach to identify features significantly associated with case-control status is conditional logistic regression (hereafter denoted as CLR) modeling. However, for high dimensional data sets, CLR can not only become computationally intensive, but can also quickly run into model convergence problems. To address these computational issues, Balasubramanian et al. [1] develop a random penalized CLR (hereafter denoted as R-PCLR) variant that merges ridge penalized CLR [9] with Random Forests [10] to identify relevant features and two-way interactions. In addition, Qian et al. [11] develop two selection methods based on the conditional and unconditional logistic likelihood functions, as well as the lasso and elastic net penalties [12, 13]. Their first method uses a two-stage approach to estimate the logistic regression parameters that are subsequently used to predict case-control status, while their second method simultaneously computes both the regression coefficient estimates and the predicted values for case-control status.

Alternately, one can approach variable selection from a Bayesian standpoint, which has several important benefits. With regards to variable selection, penalized methods identify features for inclusion by determining which of them have nonzero model coefficient estimates. Bayesian variable selection (BVS) provides more information by giving not only coefficient estimates, but also inclusion probability estimates for each feature. One common BVS approach utilizes the spike and slab prior [14, 15, 16], which involves assigning hierarchical priors to the regression coefficients by introducing indicator variables to determine whether each predictor should be considered for inclusion or removal from the model. The use of this prior is fairly widespread among BVS techniques [17, 18, 19, 20, 21], mainly due to its flexibility and ease of application, particularly for high dimensional data [22]. Lee et al. [23], Sha et al. [24], and Zhou et al. [25] all developed BVS approaches based on spike-and-slab priors [14, 15, 16] for binary outcomes, with applications to genetic microarray data. While [23] and Sha et al. [24] developed their approach using multinomial probit models by introducing latent variables, Zhou et al. [25] used logistic models. BVS methods have also yielded relatively high selection accuracy for linear regression [16, 26], and have been shown to efficiently handle ultra-high dimensional data sets [27]. Another key benefit of BVS is that it can naturally incorporate auxiliary information regarding different spatial, network, or other correlation or grouping structures among features. For instance, Smith and Fahrmeir [28] incorporate spatial correlation among features with direct applications to imaging data, while Stingo et al. [29] account for membership in a particular genetic pathway and the relationship among genes in that pathway.

To account for matching in BVS, we can specify the likelihood based on a CLR model. Based on the key benefits of BVS and CLR modeling, we propose a new methodology that formulates BVS in a CLR framework (hereafter denoted as BVS CLR) and evaluate its performance using simulation and actual studies. For comparative purposes, we also examine the performance of CLR using the lasso penalty (denoted as lasso CLR). In our applications to actual studies, we also assess the performance of BVS CLR relative to that of R-PCLR and the methods proposed by Qian et al. [11], both of which were based on lasso penalized conditional logistic likelihood functions. Several penalties exist for variable selection, e.g., elastic net, group lasso [30, 31], and sparse group lasso [32]. Since the variable selection method we propose does not take the correlation or group structure among features into account, it is more comparable to lasso. In addition, we compare our approach with that of the lasso penalty because of its ease of implementation for CLR using available software (e.g., R) relative to other variable selection penalties, e.g., SCAD [33] or adaptive lasso [34].

In Section 2, we specify the CLR model for paired data and describe our selection approach. We evaluate the performance of BVS CLR relative to lasso CLR using simulation studies in Section 3. In Section 4, we assess the performance of BVS CLR relative to lasso CLR, as well as R-PCLR and the selection approach of Qian et al. [11] which account for matching, using a matched brain imaging study of hospital acquired pneumonia (HAP) among stroke patients in Massachusetts General Hospital (MGH). We then perform a similar assessment in Section 5 using a matched study of biomarkers for near-term cardiovascular events, where we compare the performance of BVS CLR with that of lasso CLR and R-PCLR. In Section 6, we conclude with a discussion.

2 Bayesian Variable Selection for Paired Case-Control Data

Consider I case-control pairs, where Xij=(Xij1,…,XijK) denotes the observed feature values and Zij denotes case-control status for the jth member of the ith pair (i=1,…,I;j=1,2), so that Zij=1 for cases and 0 for controls. CLR models the probability that the first member of pair i is a case, given (Xi1,Xi2) and Zi1+Zi2=1, as follows:

(1)pi1=P(Zi1=1|Zi1+Zi2=1,Xi1,Xi2)=1+exp−∑k=1Kβk(Xi1,k−Xi2,k)−1,

where βk denotes the coefficient or log odds ratio for feature Xk. For features with an effect on case-control status, i.e., relevant features, βk is nonzero. From eq. (1), the conditional log-likelihood is given by

(2)lC(β)=log(LC(β))=log∏i=1Ipi1Zi1.

We now introduce γ=(γ1,…,γK), a binary vector where γk is either 1 or 0 based on whether Xk is retained in eq. (1). We assign a mixture of normal and point mass priors to the coefficient vector β=(β1,…,βK), where βk|γk∼γkN(0,σ2)+(1−γk)δ0, δ0 corresponds to π(βk=0)=1 and σ2>0. In addition, γk|ω∼Bernoulli(ω), where ω∼Beta(c,d) and c,d>0.

Based on the prior assumptions and the conditional likelihood in eq. (2), the posterior distribution of γ and β is

(3)p(β,γ|X,Z)∝LC(β)⋅π(β|γ)⋅π(γ),

such that π(β|γ)=∏k=1Kπ(βk|γk) and π(γ)=∏k=1K∫01π(γk|ω)⋅π(ω)dω.

We use Markov chain Monte Carlo (MCMC) sampling via the Metropolis-Hastings (MH) algorithm to estimate the distribution in eq. (3). Starting from random initial values for β and γ, we apply the following moves at each of S MCMC iterations:

  1. Move 1

    1. Add or remove Xr (r∈1,…,K) by choosing γr at random and changing its state, i.e., γr∗=1−γr(t), where γr∗ and γr(t) denote the proposed and current values of γr. If γr∗=1, generate the proposed βr∗ value for Xr from N(0,τ1σr∗2), where τ1(gt;0) is a proposal tuning parameter and σr∗2 is the estimated variance corresponding to the MLE of βr from a univariate CLR model on Xr. Otherwise, let βr∗=0.

      Based on the posterior probabilities p(β,γ|X,Z) and proposal densities q(⋅|β,γ) and q(⋅|γ), compute the acceptance ratio

      (4)A=p(β*,γ*|X,Z)p(β(t),γ(t)|X,Z)⋅[q(γ(t)|γ*)q(γ*|γ(t))⋅q(β(t)|β*,γ(t))q(β*|β(t),γ*)].

      Accept the proposed values (β∗,γ∗) with probability min(A,1) and retain the current values (β(t),γ(t)) otherwise.

  2. Move 2

    1. For each included Xk, i.e., Xk with γk(t)=1, generate the proposed βk∗ value from N(βk(t),τ2σk∗2), where βk(t) denotes the current value of βk, τ2(>0) is a proposal tuning parameter, and σk∗2 is the variance estimate for the univariate CLR MLE of βk.

    2. Since we do not update γ and only update the β values for the included Xk, the acceptance ratio in eq. (4) reduces to

      A=p(β∗,γ(t)|X,Z)p(β(t),γ(t)|X,Z)=LC(β∗)⋅π(β∗|γ(t))LC(β(t))⋅π(β(t)|γ(t)).

      Accept the proposed values β∗ with probability min(A,1) and retain the current values β(t) otherwise.

For each move type, we neither wanted to have proposal variances that are too small to avoid encountering mixing issues and high autocorrelations, nor did we want to have variances that are too large to avoid low acceptance rates. According to Gelman et al. [35] and Roberts et al. [36], asymptotically optimal acceptance rates for random walk Metropolis algorithms are approximately equal to 25%. Therefore, τ1 and τ2 are chosen to ensure that the acceptance rates for each move type are between 20% and 30%. Alternately, we could have used an adaptive MH algorithm, as in Lamnisos et al. [37]. However, the standard MH algorithm that we specify has good convergence properties and is easy to implement.

We then obtain the sequence (β[1],γ[1]),…,(β[S],γ[S]). Assuming a burn-in period of B iterations, estimates of the posterior inclusion probabilities p(γk|X,Z) and coefficients βk are given by

(5)pˆ(γk=1|X,Z)=∑υ=B+1Sγk[υ]S−B,βˆk=pˆ(γk=1|X,Z)⋅∑υ=B+1Sγk[υ]βk[g]∑υ=B+1Sγk[υ].

A similar averaging approach can be used to obtain variance estimates of βˆk.

Since we use standard spike-and-slab priors and Metropolis-Hastings moves to update the model parameters, the ergodicity of our MCMC sampler is guaranteed. The introduction of Move 2 does not compromise ergodicity; it is used to provide faster convergence by refining the parameter space within the selected model. To determine the convergence of β and γ in our later applications of BVS CLR to actual data, we consider a multivariate version of Gelman and Rubin’s potential scale reduction factor using the coda R package, where values substantially above 1 indicate lack of convergence [38, 39].

3 Simulation Studies

We first assess the operating performance of BVS CLR using several simulated datasets. In doing so, we focus on the level of accuracy in identifying the relevant and non-relevant features and in predicting case-control status. We also assess coefficient estimation accuracy by examining the mean squared error (MSE) of the βk estimates. In our study, we examine the performance of BVS CLR relative to that of lasso CLR with respect to selection and prediction accuracy. Although variance estimation approaches for lasso regression exist [40, 41], we take a more straightforward approach and consider, for both BVS and lasso CLR, a rough MSE approximation using the empirical mean of the squared deviations between the estimated and actual βk values. We run all analyses using R software version 3.1.2 [42].

3.1 Simulation Designs

3.1.1 Simulation of Paired Response and Feature Data

First, we simulate M=10,000 observations for variables we term as age and gender. Gender values are generated from a Bernoulli(0.5) distribution. Age values are generated from a truncated normal distribution in the range (0,60) with mean 30 and variance 100, and then rounded to the nearest integer. To examine the performance of BVS CLR for different data types, we consider the cases of binary and normal features. To simulate X values such that the first L (L<K) features are related to gender and standardized age and the first Q (Q<L<K) features are relevant, we use the following approach:

  1. Binary case Assume pm,k=P(Xk=1)=1+expage.stdm+1.5genderm−1 for observation m (m=1,…,M) for k=1,…,L and pm,k=0.5 otherwise, where age.stdm and genderm denote standardized age and gender for observation m. Each age value is standardized by subtracting its mean value and dividing by its standard deviation. Since age is included as a linear effect in our model for simulating case-control status pm,k, we standardize age to ensure that the regression coefficient we use corresponds to a reasonable effect size. To build correlation among (Xm,1,…,Xm,K), we specify correlation matrix Σ with entries ρij obtained from the phi coefficients computed from the MGH brain imaging data, where a phi coefficient measures the association between a pair of binary features. We consider the cases of both high correlation (ρij obtained from the phi coefficients ranked in the top 50 in absolute value for X1,…,XQ; ρij obtained from the remaining phi coefficient values for XQ+1,…,XK; X1,…,XQ are uncorrelated with XQ+1,…,XK) and low correlation (for X1,…,XK, ρij obtained from the phi coefficient values not ranked in the top 50 in absolute value) among the K features. Observations (Xm,1,…,Xm,K) are simulated using the bindata R package.

  2. Normal case Using the mvtnorm R package, we generate Xm,k from a normal distribution with mean age.stdm+1.5⋅genderm for k=1,…,L and 0 otherwise and covariance matrix Σ with entries ρij equal to the correlation coefficients computed from the matched biomarker study discussed in Section 5. As in the binary case, we consider the cases of both high correlation (ρij obtained from the correlation values ranked in the top 50 in absolute value for X1,…,XQ; ρij obtained from the remaining correlation values for XQ+1,…,XK; X1,…,XQ are uncorrelated with XQ+1,…,XK) and low correlation (for X1,…,XK, ρij obtained from the correlation values not ranked in the top 50 in absolute value) among the K features.

Case-control status Zm is generated from a Bernoulli(φm) distribution with φm=P(Zm=1|Xm)=1+exp−∑k=1KβkXk−0.3⋅cen.agem−0.5⋅genderm−1, where cen.agem denotes the value of age for observation m that has been centered to have mean 0. Here, β1,…,βQ fall in the range [1,2] in magnitude in the binary case, and [1,1.3] in the normal case. For each data type, values for β1,…,βQ, which capture effect sizes for the relevant features on the log odds scale, were chosen to ensure obtaining realistic odds ratios. In the normal case, BVS CLR did not converge when β1,…,βQ fell in the range 1,2. Our predictors in this case were not standardized and had a fairly wide range of values, which leads to large values for the linear predictors and, thus, a log-likelihood that diverges. Therefore, we decreased the range of values in the normal case, relative to the binary case. All remaining β elements are set to zero in both cases. From this population, we randomly select I=50 or I=200 observations with Zm=1 as cases and match them with observations with Zm=0 (controls) based on age and gender, and let these (Xij,Zij) observations constitute the training set. Omitting this set, we randomly select N=2,000(Xn,Zn)(n=1,…,N) observations from the remaining population to constitute the test set. We also consider the following scenarios:

  1. Scenario 1: K=20 features, where Q=2,5,10 are relevant

  2. Scenario 2: K=100 features, where Q=10,25,50 are relevant

  3. Scenario 3: K=600 features, where Q=60,150,300 are relevant

For each data type (binary/normal), correlation level, and scenario, we simulate 100 datasets for scenarios 1 and 2, and 20 datasets for scenario 3 due to the amount of computation time involved. Each dataset consists of a training set and test set as previously discussed. We now describe how BVS and lasso CLR are applied to each dataset.

3.1.2 Application of BVS CLR

We first apply BVS CLR to each training set, where we set the prior variance σ2 for βk|γk to 1. Since the estimated standard errors for the univariate CLR coefficient estimates were all in the 10−2 range in this study, a value of 1 for σ2 is fairly noninformative in this case. To ensure that a sufficient number of features are considered for inclusion, we assign a Beta(5,5) prior to ω in all scenarios. Although we omit their results in our discussion, we note that preliminary analyses have demonstrated the robustness of BVS CLR to specification of both σ2 (values considered ranged from 0.1 to 10) and the Beta(c,d) prior.

We run BVS CLR for S=50,000 iterations (B=20,000 burn-in) in scenario 1, S=80,000 iterations (B=30,000 burn-in) in scenario 2, and S=100,000 iterations (B=50,000 burn-in) in scenario 3. Along with the estimates pˆ(γk=1|X,Z) and βˆk in (5), we compute the variance estimates of βˆk. For each post burn-in iteration υ (υ=B+1,…,S), we compute the estimated case probability pn[υ]=1+exp−∑k=1Kβk[υ]Xn,k−1 for the nth test set observation. Using pn[υ], we compute the Bayesian model averaged (BMA) case probability pˆn,BMA using the following approaches: (1) average pn[B+1],…,pn[S], (2) average Zn[B+1],…,Zn[S], where Znυ=1 if pn[υ]≥0.5 and 0 otherwise, and (3) generate Znυ from a Bernoulli(pn[υ]) distribution and average Zn[B+1],…,Zn[S]. We discuss the results obtained using the first approach, although all three approaches give similar results.

3.1.3 Application of Lasso CLR

We also examine the performance of lasso CLR, which is based on the penalized conditional log-likelihood

(6)lC(β)−λ∑k=1K|βk|.

In eq. (6), the feature values are standardized and λ≥0 is a tuning parameter that is typically estimated using V-fold cross validation. Use of the lasso penalty in eq. (6) shrinks all βk estimates to zero, yielding a subset of features found to be relevant in eq. (1) due to having nonzero βk estimates. Using the survival and penalized R packages, we apply 10-fold cross-validation to estimate λ and run lasso CLR on each training set, yielding coefficient estimates βˆk,las and case probability estimates pˆn,las=1+exp−∑k=1Kβˆk,lasXn,k−1 for the nth test set observation.

A schematic diagram summarizing the simulation details, and BVS and lasso CLR applications, is provided in Figure 1.

Figure 1 Schematic diagram of simulation and application details for each simulation.
Figure 1

Schematic diagram of simulation and application details for each simulation.

3.2 Variable Selection Accuracy

To measure selection accuracy for BVS and lasso CLR, we compute the areas under the ROC curves (AUCs) using inclusion probability estimates pˆ(γk=1|X,Z) for BVS CLR and the magnitudes of the coefficient estimates |βˆk,las| for lasso CLR. The AUC was used as a metric for selection accuracy because it gives an overall picture of the degree to which we can correctly identify relevant and non-relevant features. However, we acknowledge that, in practice, we would have to rely on the use of a defined threshold for identifying relevant and non-relevant features.

We present the median and interquartile range (IQR) values for these AUCs across simulations in the normal and binary cases for scenarios 1 and 2 in Appendix Table 5 and for scenario 3 in Appendix Table 6. For both approaches, we have that, regardless of the number of pairs I, data type, and correlation level among the relevant features, selection accuracy decreases as the total number of features K increases for a given number of relevant features Q, and as Q increases for a given K. We also have that for both approaches, selection accuracy increases with increasing I for a given Q and K, regardless of data type and correlation level. Also, as the level of correlation among the relevant features increases from low to high, selection accuracy remains comparable (magnitude of difference less than 0.05) or decreases in both the normal and binary cases. In cases where selection accuracy decreases, we have that the decrease is generally more pronounced for lasso CLR.

Regardless of the percentage of relevant features, BVS CLR generally yields higher selection accuracy than lasso CLR in the normal and binary cases when K=100 features and I=50 pairs. We also note an improved performance of BVS CLR for normally distributed features with high correlation when K=100 features and I=200 pairs. With K=20 features, BVS CLR yields higher accuracy for the binary high correlation case with I=50 regardless of the percentage of relevant features. When 25% and 50% of the features are relevant, this improved accuracy is also observed for the binary low correlation case and for the normal high correlation setting. With K=600 features, the two methods have a comparable performance.

3.3 Prediction Accuracy

We also examine prediction accuracy for BVS and lasso CLR, which we measure using the AUCs based on the predicted case probabilities pˆn,BMA and pˆn,las. For each approach, the median and IQR values for these AUCs across simulations are presented in Appendix Table 8 for scenarios 1 and 2 and in Appendix Table 9 for scenario 3.

In general, prediction accuracy is higher for the normal datasets than for the binary datasets, and also increases as the number of pairs I increases. Also, prediction accuracy is highest when 25% of the features are relevant for K=20,100, and remains comparable across the number of relevant features Q for K=600.

BVS CLR has higher prediction accuracy than lasso CLR for the simulation setting with I=50 pairs when K=100 features and 50% of them are relevant, for both normal and binary cases with low or high correlation. We also observe improved prediction accuracy of BVS CLR for K=600 features with 50% relevant in all cases (normal, binary, low or high correlation, I=50 or 200 pairs). For K=600 with 25% relevant, BVS CLR has higher prediction accuracy in all cases when I=50 pairs, and for the normal setting when I=200 pairs. For K=600 features with 10% relevant, BVS CLR has higher prediction accuracy for the normal low correlation case with I=50 pairs. For all other settings, we find comparable prediction performance between BVS CLR and lasso CLR.

3.4 MSE for Coefficient Estimates

For each scenario, we then assess the level of coefficient estimation accuracy for BVS and lasso CLR by examining the average MSE for each feature across simulations, and then computing the median of these averaged MSEs across the relevant features and the non-relevant features. These results are reported in Appendix Table 11 and Table 13 for the relevant and non-relevant features, respectively.

Relative to lasso CLR, these median MSEs are generally lower for BVS CLR with respect to the relevant features and comparable with respect to the non-relevant features regardless of data type, correlation level among the relevant features, number of pairs, number of features (total and relevant). With fewer pairs, this decrease in MSE for BVS CLR with respect to the relevant features is more pronounced. For both approaches, we see that MSEs for the relevant features are generally larger for binary datasets, increase as the total number of features K and number of relevant features Q increase, and decrease as the number of pairs I increases. However, no apparent pattern emerges in our results as we increase the correlation level among the relevant features.

Although we do not present the results in our discussion, the level of estimation accuracy for BVS CLR was also assessed using coverage probabilities based on the 95% highest posterior density (HPD) intervals for each coefficient. In doing so, we observed high coverage probabilities across features and HPD intervals with widths that grew narrower with increasing I.

3.5 Convergence of BVS CLR

To explore the convergence performance of BVS CLR, we examine, across post burn-in iterations, trace plots of the log posterior probabilities in Figure 5 and the number of selected features in Figure 7 for the normal case when I=50 pairs and 25% of the features are relevant and weakly correlated. Based on these plots, we do not see any evidence of non-convergence. Although we omit their results, similar patterns are found in the trace plots for all other simulation scenarios, data types, number of pairs, and levels of correlation among the relevant features. To explore the behavior of features with the highest inclusion probabilities, we consider, as an example, the case of K=20 normal features of which Q=5 are relevant and weakly correlated. For this case, we examine in Figure 8 a trace plot of the inclusion status (yes/no) across post burn-in iterations in one simulation for features X1,X3,X4 whose inclusion probabilities are ranked among the top three for I=50 pairs.

3.6 Accuracy under Reduced Coefficient Values

To assess the performance of BVS CLR relative to lasso CLR in the case when β1,…,βQ for the relevant features are relatively small, we also consider the case where β1,…,βQ fall in the range [0.3,0.7] in magnitude, in scenarios 1 (K=20 features) and 2 (K=100 features). For both the binary and normal cases, the median and IQR values for the selection and prediction accuracy AUCs across simulations are reported in Appendix Table 7 and Table 10, respectively. In Appendix Table 12 and Table 14, we report the medians of the averaged MSEs across the relevant features and across the non-relevant features, respectively. In summary, relative to when β1,…,βQ are at least 1, the improvement in performance for BVS CLR compared with lasso CLR is more pronounced with respect to selection accuracy, relatively unchanged with respect to prediction accuracy and MSE for the non-relevant features, and less pronounced with respect to MSE for the relevant features.

4 Application to MGH HAP Imaging Study

4.1 Description

We first examine a case-control brain imaging study conducted by Kemmling et al. [43], in which acute ischemic stroke patients admitted to the Stroke Service Unit at MGH stroke service were classified according to whether or not they met the criteria for having hospital acquired pneumonia (HAP), i.e., suspicion or mention of pneumonia in the patient’s medical record at least 48 hours after admission requiring antibiotic treatment. 215 acute ischemic stroke patients classified as having HAP were then matched with 215 non-HAP acute ischemic stroke patients on the basis of age, gender, and NIH stroke scale (NIHSS) upon admission. In this study, each patient was measured on both clinical and neuroimaging features. Clinical features include age, gender, admission NIHSS, length of hospitalization, as well as the presence/absence of the following: dysphagia, dyslipidemia, smoking history, coronary artery disease, diabetes mellitus, atrial fibrillation, hypertension, and in-hospital mortality. In Table 1, we present the summary statistics for each of these features.

Table 1

Clinical features of matched HAP and non-HAP patients [mean ± standard error for age and length of hospitalization; median (interquartile range) for admission NIHSS; n (%) for the remaining features].

Clinical FeatureHAP (n=215)non-HAP (n=215)p-Value
Age (years)a72.2±14.972.3±13.9–
Malea116 (54%)116 (54%)–
Admission NIHSSa13 (6−19)13 (6−19)–
Dysphagia151 (70.2%)151 (70.2%)1.00b
Hypertension146 (68.0%)139 (64.7%)0.53b
Dyslipidemia74 (34.4%)68 (31.6%)0.59b
Diabetes mellitus51 (23.7%)44 (20.5%)0.49b
Atrial fibrillation65 (30.2%)58 (27.0%)0.51b
Smoking history38 (17.7%)35 (16.3%)0.78b
Coronary artery disease59 (27.4%)51 (23.7%)0.39b
Mortality41 (19.1%)38 (17.7%)0.80b
Length of hospitalization (days)12.8±10.26.1±4.6<0.0001c

  1. Feature used to match HAP and non-HAP patients.

  2. McNemar’s test used to compare HAP and non-HAP patients.

  3. Wilcoxon signed-rank test used to compare HAP and non-HAP patients.

To extract the neuroimaging features for each patient, subacute ischemic brain lesions were first outlined slice-by-slice in diffusion weighted magnetic resonance imaging (MRI-DWI) or computerized tomography (CT) images, each of which were chosen with an acquisition time approximately 48 hours after symptom onset. MRI-DWI/CT images and their respective binary lesion masks were affine registered to standard MNI-152 space and manually corrected for registration errors. The same imaging protocol was used for all patients.

All lesion masks were then segmented into 68 pairs (left-right hemispheres) of cortical and subcortical/brainstem white matter brain regions based on the “Johns Hopkins University white-matter” and “Harvard-Oxford cortical structural” atlases, which were created by standardized anatomic labeling of multiple subjects linearly registered to MNI-152 standard space [44, 45]. Binarized atlases defining a specific structure with at least 25% probability of anatomic localization were used. For all patients, the percentage of infarction in a specific brain region was first measured and then dichotomized using its median value as having zero or positive infarction. Specifically, the neuroimaging feature examined for each brain region was the presence/absence of positive infarction in that region. Given the relatively small sample size, brain regions with positive infarction for less than 5% of HAP patients or non-HAP patients will yield unstable log odds ratio estimates [11]. Therefore, to stabilize calculations, we only consider the 130 brain regions that have positive infarction for at least 5% of HAP patients and 5% of non-HAP patients. Another examined neuroimaging feature was infarction volume. To avoid model fitting issues due to its highly skewed distribution and the presence of outliers, we categorized infarction volume using the tertiles of its distribution. Two indicator variables were introduced to indicate whether infarction volume was at least equal to its first tertile, and whether it was at least equal to its second tertile.

In applying our proposed methodology, we aim to identify the subset of clinical and neuroimaging features that are most associated with, or most relevant to, having HAP. Although there is not necessarily a link between HAP and white matter infarctions, per se, associations between HAP and infarctions in specific brain regions have been found in prior studies. Kemmling et al. [43] showed that infarction in the right hemispheric peri-insular cortical regions was associated with the risk of acquiring HAP in acute ischemic stroke patients. They explain this finding by discussing the fact that previous studies have demonstrated the association of the right insular region with autonomically-induced immunosuppression and susceptibility to infection [46, 47, 48], as well as the association of right hemispheric peri-insular infarction to autonomic dysfunction and pathologic sympathetic activity [49].

In this application, we evaluate the performance of BVS CLR, and also assess its performance relative to lasso CLR, R-PCLR, and the selection approach of Qian et al. [11].

4.2 Method

We run BVS and lasso CLR in 50 parallel chains, i.e., both methods are applied to the data set 50 times, where a different random seed is used in each instance. To each chain, we apply (1) BVS CLR to obtain inclusion probability estimates pˆ(γk=1|X,Z), and (2) lasso CLR to obtain coefficient estimates βˆk,las, where 10-fold cross-validation is used to estimate the tuning parameter λ in eq. (6). In applying BVS CLR, we use S=80,000 iterations (B=40,000 burn-in) for each chain. As in our simulation study, we set the variance σ2 in our normal prior for βk|γk to 1 and assign a Beta(5,5) prior to ω. We retained these values due to the robustness of BVS CLR to the specification of σ2 and the Beta(c,d) prior that we observed in preliminary analyses for our simulation study.

For lasso CLR, we identified features with nonzero coefficient estimates as relevant. To determine the threshold value for identifying relevant features in BVS CLR, we used the posterior inclusion probability estimates pˆ(γk=1|X,Z) and applied the Bayesian FDR approach proposed by Muller et al. [50] and used in prior studies [51, 52]. This is a variation of the Benjamini and Hochberg [53] procedure where the threshold is based on increments in the ordered posterior probabilities rather than ordered p-values. We used an FDR level of 0.27 in order to identify approximately as many as features as in lasso CLR.

4.3 Results

In Table 2 and Table 3, we present the following, averaged across the 50 parallel chains: (1) BVS CLR inclusion probabilities for the features identified using the Bayesian FDR approach and their corresponding ranks, along with their coefficient estimates and standard deviations (SDs) and, (2) the lasso CLR coefficient estimates and their SDs for the features identified as relevant, along with their ranks based on the magnitude of their coefficient estimates.

There is 64% overlap among the features identified as relevant under BVS and lasso CLR. There is also a considerable degree of selection overlap using the approach of Qian et al. [11], which identifies as relevant 71% of the features selected by BVS CLR and 85% of those selected by lasso CLR. On the other hand, only 21% of the features selected by BVS CLR and 36% of the features selected by lasso CLR are also identified as relevant using R-PCLR. Relative to lasso CLR, the coefficient estimates under BVS CLR are noticeably larger in magnitude. Although the SDs for the coefficient estimates are smaller under lasso CLR, this may be due to the fact that the lasso CLR coefficient estimates are generally close to zero in magnitude. In addition, the multivariate scale reduction factors for β and γ were 1.15 and 1.12, respectively. Since these values are not substantially greater than 1, these results do not indicate a lack of convergence.

Table 2

(MGH Imaging Study) BVS inclusion rates (and corresponding ranks), coefficient estimates and SDs for BVS CLR.

CodeRegion/VolumeRankInclusion Prob.Coef. Est.Coef. SD
v67a,bInfarction volume ≥67th percentile10.860.900.13
j10aCerebral peduncle R20.850.980.10
j34aFornix (cres) / Stria terminalis R30.820.970.13
h74aTemporal Fusiform Cortex – anterior division R40.82−1.120.14
h46bLateral Occipital Cortex – inferior division R60.740.770.12
j23bPosterior thalamic radiation L80.67−0.570.13
h32aInferior Temporal Gyrus – temporooccipital part R100.670.630.11
h51aJuxtapositional Lobule Cortex L110.66−0.640.12
h5aSuperior Frontal Gyrus L130.64−0.540.11
j19aSuperior corona radiata L50.760.700.15
h72aLingual Gyrus R70.730.700.10
cadaCoronary artery disease90.670.450.06
j30Cingulum (cingulate gyrus) R120.65−0.550.09
h33Postcentral Gyrus L140.640.520.13

  1. Selected using Qian et al. [11] selection approach

  2. Selected using R-PCLR.

Table 3

(MGH Imaging study) Cross-validated coefficient estimates and SDs (and corresponding ranks for magnitude of coefficient estimates) for lasso CLR.

CodeRegion/VolumeRankCoef. Est.Coef. SD
v67a,bInfarction volume ≥67th percentile40.180.02
j10aCerebral peduncle R60.120.02
j34aFornix (cres) / Stria terminalis R20.270.01
h74aTemporal Fusiform Cortex – anterior division R10−0.010.02
h46bLateral Occipital Cortex – inferior division R11<0.010.01
j23bPosterior thalamic radiation L14<−0.01<0.01
h32aInferior Temporal Gyrus – temporooccipital part R10.420.01
h51aJuxtapositional Lobule Cortex L3−0.260.02
h5aSuperior Frontal Gyrus L8−0.07<0.01
h8a,bMiddle Frontal Gyrus R50.130.01
j26a,bSagittal stratum R70.08<0.01
v33aInfarction volume ≥33rd percentile90.010.02
h45aLateral Occipital Cortex – inferior division L12<-0.01<0.01
h25aMiddle Temporal Gyrus – temporooccipital part L13<0.01<0.01

  1. Selected using Qian et al. [11] selection approach.

  2. Selected using R-PCLR.

Figure 2 and 4.4 display heatmaps of Cramèr’s V matrix plot for the pairwise differences among the features selected in Table 2 and Table 3 as relevant and those not selected, where the pairwise feature differences are considered to account for matching. In these plots, black denotes perfect positive/negative correlation, white denotes no correlation, and light/dark gray denotes weak/strong correlation. In both plots, we see that the correlation among features is generally low, so that the structure of this dataset most closely corresponds to Scenario 2 (K=100 features) in our simulation study when all binary features were weakly correlated.

Most of the right hemispheric brain regions identified under BVS CLR were found to be associated with HAP in Kemmling et al. [43] , who showed that infarction in the right hemispheric peri-insular cortical regions was associated with the risk of acquiring HAP in acute ischemic stroke patients. As discussed in Section 4.1, one explanation Kemmling et al.provided for this finding was the association of the right insular region with autonomically-induced immunosuppression and susceptibility to infection shown in Cechetto and Chen [46], Sander and Klingelhöfer [47], and Meyer et al. [48], as well as the association of right hemispheric peri-insular infarction to autonomic dysfunction and pathologic sympathetic activity shown in Colivicchi et al. [49].

Figure 2: (MGH Imaging Study) Cramer’s V matrix plots of paired differences for features selected (a) and not selected (b) in Table 2 and Table 3 (black denotes perfect positive or negative correlation, white denotes no correlation; light/dark gray denotes weak/strong correlation).
Figure 2:

(MGH Imaging Study) Cramer’s V matrix plots of paired differences for features selected (a) and not selected (b) in Table 2 and Table 3 (black denotes perfect positive or negative correlation, white denotes no correlation; light/dark gray denotes weak/strong correlation).

5 Application to Cardiovascular Disease Biomarker Study

5.1 Description

The cardiovascular disease biomarker study we examine was a matched case-control study conducted by the High Risk Plaque Initiative [BG Medicine Inc.(Waltham, MA) and other partners] to discover prognostic biomarkers in blood plasma for near-term cardiovascular events. Subjects from the CATHGEN study were selected for this investigation. The CATHGEN project collected peripheral blood samples from consenting research subjects undergoing cardiac catheterization at Duke University Medical Center from 2001 to 2011. 68 cases were selected from among individuals who had a major adverse cardiac event (MACE) within two years following the time of their sample collection. In a 1:1 matched study design, 68 controls were selected from individuals who were MACE-free for the two years following sample collection and were matched to cases on age, gender, race/ethnicity and severity of coronary artery disease. High-content mass spectrometry and multiplexed immunoassay-based techniques were employed to quantify 625 proteins and metabolites from each subject’s serum specimen. Comprehensive metabolite profiling of the individual samples was based on a combination of four platforms employing mass spectrometry (MS) based techniques to profile lipids, fatty acids, amino acids, sugars and other metabolites. Proteomic analysis was based on a combination of targeted methods using a quantitative multiplexed immunoassay technique as well as a comprehensive protein profiling strategy based on tandem mass spectrometry. A detailed description of the mass spectrometry based platforms and proteomics analysis can be found in a previous publication [54]. In our analysis, the identities of the measured metabolites and proteins are masked due to a data confidentiality agreement with the stakeholders involved in the study.

We consider for analysis the 593 biomarkers with complete data. In this application, we evaluate BVS CLR under the assumption of normality for all biomarkers, and assess its performance in relation to lasso CLR and R-PCLR.

5.2 Method

We used the same procedure described in Section 4.2, except that 25 parallel chains with different random seeds are run. To account for the high dimensionality of this dataset, we use S=150,000 iterations (B=75,000 burn-in) for each chain. As in our previous application, we set the variance σ2 in our normal prior for βk|γk to 1 and assign a Beta(5,5) prior to ω. To identify relevant features, we use the Bayesian FDR approach (at an FDR level of 0.41) for BVS CLR based on the inclusion probability estimates and examine features whose lasso CLR coefficient estimates are nonzero.

5.3 Results

In Table 4, we present the estimates of the BVS CLR inclusion probabilities across the 25 chains for the biomarkers identified using the Bayesian FDR approach and their corresponding ranks, along with the corresponding coefficient estimates and their SDs. We also report the lasso CLR coefficient estimates and their SDs across the chains for the biomarkers identified as relevant under lasso CLR, along with their ranks based on the magnitude of their coefficient estimates.

There is 52% overlap among the biomarkers identified as relevant under BVS and lasso CLR, and only 25% overlap among the biomarkers identified as relevant under BVS CLR and R-PCLR. Relative to lasso CLR, the coefficient estimates under BVS CLR are larger in magnitude, along with their SDs, which may be due to the fact that the lasso CLR coefficient estimates are generally close to zero.

In Figure 3 and Figure 4, we display heatmaps of the correlation matrix plot for the pairwise differences among the biomarkers selected in Table 4 as relevant and those not selected, where we see that the correlation among biomarkers is generally low, so that the structure of this dataset most closely resembles Scenario 3 (K=600 features) in our simulation study when all normal features were weakly correlated. The fact that we saw lower MSE in our application of BVS CLR to these normal datasets suggests that the coefficient estimates under BVS CLR are more likely to be more accurate relative to lasso CLR.

We note that the multivariate scale reduction factors for β and γ were 3.01 and 2.95, respectively. Although we do not report the results, increasing the number of iterations and folds to a reasonable degree, considering the computational expense involved, did not substantially decrease these reduction factor values. This may result from the high number of biomarkers examined, relative to the number of subjects.

Table 4

(CVS Disease Biomarker Study) BVS CLR inclusion rates and corresponding ranks; BVS and lasso CLR CV coefficient estimates (and SDs).

BVS CLR Lasso CLR
BiomarkerRankInclusion Prob.Coef. Est.Coef. SDBiomarkerRankCoef. Est.Coef. SD
V59510.68−0.760.15V5953−0.270.05
V27220.670.780.18V27250.180.06
V58230.64−0.660.15V5829−0.130.01
V275a40.630.640.19V275a120.100.01
V617a50.63−0.630.28V617a2−0.410.05
V16460.62−0.580.13V1644−0.190.02
V625a70.610.630.18V625a10.540.06
V17480.610.520.17V17460.170.03
V59390.610.560.22V593100.110.03
V436100.60−0.510.25V43626−0.010.01
V362110.590.510.16V362130.090.04
V607120.580.460.14V60770.160.04
V535150.580.450.18V535180.060.03
V75230.57−0.440.20V7515−0.080.03
V31260.560.490.29V3128<0.010.01
V464300.56−0.370.22V4648−0.150.05
V158130.58−0.350.20V185110.100.01
V219140.580.430.17V66140.080.01
V620160.58−0.340.18V314160.080.03
V24170.580.460.17V122170.080.01
V151180.57−0.450.20V605a190.030.01
V96190.570.410.13V41200.020.01
V223200.570.220.23V611a210.020.02
V234210.570.280.18V563220.010.01
V227220.570.300.17V495230.010.01
V152240.56−0.380.17V45924−0.010.02
V99250.56−0.440.20V64250.010.01
V54270.560.470.14V621270.010.01
V89280.56−0.340.13V12829<0.010.01
V228290.560.360.14V11730<0.010.01
V578310.560.380.13V597a31<0.010.01

  1. Note:Selected using R-PCLR.

Figure 3: (CVS Disease Biomarker Study) Correlation matrix plots of paired differences for biomarkers selected (a) and not selected (b) in Table 4 (black denotes perfect positive or negative correlation, white denotes no correlation; light/dark gray denotes weak/strong correlation).
Figure 3:

(CVS Disease Biomarker Study) Correlation matrix plots of paired differences for biomarkers selected (a) and not selected (b) in Table 4 (black denotes perfect positive or negative correlation, white denotes no correlation; light/dark gray denotes weak/strong correlation).

6 Discussion

In examining the simulation study results for BVS CLR in Section 3, we see that selection and prediction accuracy increase as the number of pairs I increases, while MSE (for the relevant features) decreases. Both selection accuracy and MSE increase with decreasing K and Q, while selection accuracy also decreases as the correlation between features increases. Prediction accuracy is higher and MSE is lower for normal datasets, while prediction accuracy is higher for datasets where the numbers of total and relevant features are both moderate relative to the other simulation scenarios considered, i.e., K=100 where 25% of features are relevant. Relative to lasso CLR, BVS CLR generally has lower MSE. In addition, BVS CLR most often has higher selection accuracy when the total number of features is moderate (K=100) and higher prediction accuracy when we have at least a moderate number of features, which is more pronounced in datasets with fewer pairs. We note that the increase in selection accuracy of BVS CLR, relative to lasso CLR, is especially pronounced when the β values for the relevant features is relatively small (i.e., less than 1). Based on the study results in Section 4 and Section 5, the degree of overlap in features identified as relevant is higher for the binary MGH imaging data.

In cases where selection accuracy is solely of interest, we can modify BVS CLR using the data augmentation approach in Sha et al. [24] based on the probit approximation to the logit link [1+exp(−ν)]−1≈Φ(ν/1.7) [57]. This formulation allows for the integration of β out of the likelihood function, so that only the posterior inclusion probabilities are estimated. When applied to high dimensional datasets, this approximation approach can increase computational efficiency and improve mixing of the MCMC sampler [24]. Although we do not report its results, we also consider the probit approximation pi1=Φ11.7∑k=1Kβk(Xi1,k−Xi2,k) to the case probability in eq. (1) when applying our proposed BVS approach and obtain nearly identical results compared with BVS CLR.

Possible future research directions include an extension of BVS CLR to account for interactions among the examined features. This can be done for two-way interactions by modifying our selection approach based on the methodology developed by Chipman [55]. Extensions of BVS CLR that incorporate specific correlation or grouping structures among different features, including extensions of the approaches of Smith and Fahrmeir [28] and Stingo et al. [29] to incorporate the spatial correlation and the correlation among features known to be involved in similar biological functions, can be considered for matched case-control data.

Although our selection methodology focuses on 1:1 case-control matching, we can extend our formulated models for BVS CLR to handle the more general case of 1:n case-control matching. Another extension would account for matching across multiple groups when an intrinsic ordering exists among the case and control groups, e.g., matching individuals across different disease severity levels. We will further investigate how to do so using the conditional adjacent categories logistic modeling approach developed by Mukherjee et al. [56] for matched 1:n case-control studies. Through appropriate transformations of the features and notational expansion of the matched sets, we can re-frame the modeling approach of Mukherjee et al. [56] for BVS CLR to handle ordinal-based matching.

Funding statement: This work was supported by the National Institutes of Health (Grant / Award Number: T32NS048005; Grant / Award Number: F32NS081904). The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health.

Appendix

A Simulation Study Tables

Table 5

(Simulation Study: Selection Accuracy) Summary measures for AUCs obtained from BVS inclusion probabilities and magnitudes of lasso CLR coefficient estimates for K=20 and K=100 features (β1,…,βQ at least 1 in magnitude).

Scenario 1Scenario 2
20 features100 features
Median(IQR)Median(IQR)
Data typeCorrelation levelNumber of pairsBVS CLRLasso CLRBVS CLRLasso CLR
NormalaLow501(<0.01)1(<0.01)0.89(0.1)0.7(0.19)
2001(<0.01)1(<0.01)1(<0.01)1(<0.01)
High501(<0.01)1(0.08)0.77(0.12)0.61(0.07)
2001(<0.01)1(<0.01)0.87(0.07)0.68(0.07)
BinaryaLow500.94(0.12)0.92(0.41)0.81(0.1)0.69(0.13)
2001(<0.01)1(<0.01)0.99(0.02)0.99(0.01)
High500.94(0.14)0.72(0.5)0.78(0.11)0.68(0.11)
2001(<0.01)1(<0.01)0.98(0.02)0.99(0.03)
NormalbLow501(0.04)0.99(0.04)0.75(0.08)0.58(0.12)
2001(<0.01)1(<0.01)0.98(0.02)0.98(0.03)
High500.95(0.07)0.73(0.07)0.69(0.08)0.54(0.05)
2001(0.03)0.99(0.14)0.84(0.06)0.57(0.04)
BinarybLow500.91(0.12)0.87(0.16)0.73(0.07)0.64(0.1)
2001(<0.01)1(<0.01)0.96(0.03)0.95(0.04)
High500.92(0.11)0.86(0.16)0.75(0.08)0.63(0.07)
2001(<0.01)1(<0.01)0.91(0.04)0.89(0.05)
NormalcLow500.93(0.08)0.88(0.12)0.64(0.06)0.5(0.01)
2001(<0.01)1(<0.01)0.9(0.04)0.86(0.05)
High500.96(0.07)0.68(0.09)0.67(0.08)0.51(0.05)
2000.98(0.05)0.68(0.07)0.84(0.07)0.58(0.04)
BinarycLow500.87(0.13)0.78(0.16)0.66(0.07)0.55(0.06)
2001(0.02)1(0.01)0.89(0.04)0.84(0.05)
High500.84(0.12)0.74(0.14)0.67(0.08)0.55(0.05)
2000.99(0.03)0.99(0.03)0.8(0.07)0.77(0.06)

  1. 10% relevant (Q=2 for K=20; Q=10 for K=100);

  2. 25% relevant (Q=5 for K=20; Q=25 for K=100);

  3. 50% relevant (Q=10 for K=20; Q=50 for K=100).

Table 6

(Simulation Study: Selection Accuracy) Summary measures for AUCs obtained from BVS inclusion probabilities and magnitudes of lasso CLR coefficient estimates for K=600 features (β1,…,βQ at least 1 in magnitude).

Scenario 3
600 features
Median(IQR)
Data typeCorrelation levelNumber of pairsBVS CLRLasso CLR
NormalaLow500.52(0.05)0.5(0.01)
2000.61(0.05)0.65(0.06)
High500.57(0.05)0.53(0.02)
2000.63(0.04)0.57(0.02)
BinaryaLow500.53(0.03)0.51(0.02)
2000.69(0.04)0.71(0.03)
High500.54(0.06)0.52(0.02)
2000.66(0.04)0.67(0.03)
NormalbLow500.5(0.03)0.5(0.01)
2000.57(0.03)0.56(0.03)
High500.53(0.06)0.52(0.02)
2000.59(0.03)0.57(0.01)
BinarybLow500.51(0.03)0.5(0.01)
2000.61(0.05)0.58(0.02)
High500.52(0.05)0.51(0.01)
2000.62(0.04)0.6(0.03)
NormalcLow500.51(0.04)0.5(0.01)
2000.54(0.03)0.51(0.02)
High500.52(0.04)0.51(0.01)
2000.55(0.03)0.55(0.01)
BinarycLow500.51(0.03)0.5(<0.01)
2000.57(0.03)0.54(0.01)
High500.54(0.04)0.53(0.01)
2000.61(0.02)0.6(0.01)

  1. 10% relevant (Q=60).

  2. 25% relevant (Q=150).

  3. 50% relevant (Q=300).

Table 7

(Simulation Study: Selection Accuracy) Summary measures for AUCs obtained from BVS inclusion probabilities and magnitudes of lasso CLR coefficient estimates for K=20 and K=100 features (β1,…,βQ between 0.3 and 0.7 in magnitude).

Scenario 1Scenario 2
20 features100 features
Median(IQR)Median(IQR)
Data typeCorrelation levelNumber of pairsBVS CLRLasso CLRBVS CLRLasso CLR
NormalaLow500.83(0.25)0.5(0.21)0.75(0.09)0.53(0.13)
2001(<0.01)1(<0.01)0.94(0.08)0.94(0.07)
High500.89(0.22)0.5(0.05)0.68(0.14)0.53(0.08)
2001(<0.01)1(<0.01)0.8(0.1)0.66(0.08)
BinaryaLow500.64(0.28)0.5(<0.01)0.62(0.11)0.5(0.04)
2000.88(0.19)0.73(0.3)0.8(0.08)0.73(0.14)
High500.69(0.28)0.5(<0.01)0.61(0.12)0.5(0.05)
2000.83(0.21)0.71(0.28)0.82(0.08)0.74(0.12)
NormalbLow500.88(0.14)0.85(0.17)0.7(0.07)0.51(0.07)
2001(<0.01)1(<0.01)0.9(0.05)0.89(0.07)
High500.83(0.12)0.6(0.12)0.71(0.07)0.55(0.04)
2000.93(0.08)0.68(0.11)0.83(0.05)0.62(0.05)
BinarybLow500.67(0.21)0.5(0.03)0.59(0.09)0.51(0.04)
2000.84(0.11)0.68(0.17)0.78(0.09)0.74(0.07)
High500.7(0.16)0.5(0.05)0.63(0.09)0.52(0.07)
2000.89(0.09)0.83(0.14)0.75(0.1)0.71(0.07)
NormalcLow500.82(0.13)0.77(0.16)0.63(0.07)0.5(0.01)
2000.99(0.03)0.99(0.03)0.84(0.05)0.81(0.05)
High500.85(0.12)0.52(0.1)0.63(0.08)0.49(0.03)
2000.96(0.08)0.67(0.08)0.72(0.07)0.49(0.05)
BinarycLow500.71(0.15)0.5(0.11)0.61(0.07)0.51(0.04)
2000.88(0.1)0.82(0.16)0.76(0.06)0.71(0.05)
High500.7(0.13)0.5(0.07)0.6(0.06)0.51(0.04)
2000.88(0.08)0.82(0.11)0.64(0.12)0.62(0.05)

  1. 10% relevant (Q=2 for K=20; Q=10 for K=100).

  2. 25% relevant (Q=5 for K=20; Q=25 for K=100).

  3. 50% relevant (Q=10 for K=20; Q=50 for K=100).

Table 8

(Simulation Study: Prediction Accuracy) Summary measures for AUCs obtained from true case-control status and BVS CLR/lasso CLR predicted case probabilities in independent test sets for K=20 and K=100 features (β1,…,βQ at least 1 in magnitude).

Scenario 1Scenario 2
20 features100 features
Median(IQR)Median(IQR)
Data typeCorrelation levelNumber of pairsBVS CLRLasso CLRBVS CLRLasso CLR
NormalaLow500.66(0.13)0.65(0.15)0.71(0.13)0.65(0.24)
2000.67(0.08)0.66(0.09)0.79(0.07)0.79(0.1)
High500.6(0.14)0.55(0.15)0.61(0.14)0.7(0.13)
2000.64(0.08)0.63(0.09)0.72(0.08)0.73(0.07)
BinaryaLow500.55(0.09)0.51(0.07)0.6(0.1)0.61(0.16)
2000.53(0.06)0.53(0.06)0.69(0.08)0.69(0.09)
High500.53(0.11)0.5(0.05)0.63(0.08)0.66(0.13)
2000.51(0.06)0.51(0.05)0.69(0.07)0.69(0.07)
NormalbLow500.87(0.13)0.89(0.1)0.79(0.16)0.64(0.29)
2000.88(0.04)0.88(0.05)0.89(0.11)0.91(0.08)
High500.89(0.04)0.9(0.03)0.91(0.02)0.91(0.02)
2000.9(0.02)0.9(0.02)0.94(0.01)0.94(0.01)
BinarybLow500.72(0.11)0.72(0.1)0.72(0.09)0.69(0.15)
2000.73(0.06)0.74(0.06)0.83(0.05)0.84(0.06)
High500.74(0.1)0.74(0.1)0.68(0.1)0.63(0.12)
2000.75(0.03)0.75(0.04)0.8(0.05)0.79(0.06)
NormalcLow500.77(0.16)0.76(0.19)0.76(0.09)0.5(0.06)
2000.82(0.08)0.83(0.09)0.87(0.11)0.86(0.1)
High500.79(0.08)0.78(0.12)0.87(0.03)0.79(0.04)
2000.83(0.05)0.84(0.05)0.91(0.02)0.86(0.03)
BinarycLow500.69(0.1)0.67(0.17)0.69(0.07)0.58(0.12)
2000.71(0.06)0.71(0.07)0.78(0.06)0.77(0.08)
High500.63(0.1)0.61(0.14)0.66(0.09)0.58(0.1)
2000.63(0.06)0.62(0.06)0.78(0.05)0.78(0.06)

  1. 10% relevant (Q=2 for K=20; Q=10 for K=100).

  2. 25% relevant (Q=5 for K=20; Q=25 for K=100);

  3. 50% relevant (Q=10 for K=20; Q=50 for K=100).

Table 9

(Simulation Study: Prediction Accuracy) Summary measures for AUCs obtained from true case-control status and BVS CLR/lasso CLR predicted case probabilities in independent test sets for K=600 features (β1,…,βQ at least 1 in magnitude).

Scenario 3
600 features
Median(IQR)
Data typeCorrelation levelNumber of pairsBVS CLRLasso CLR
NormalaLow500.66(0.03)0.55(0.07)
2000.79(0.02)0.76(0.02)
High500.77(0.07)0.75(0.03)
2000.84(0.04)0.81(0.03)
BinaryaLow500.56(0.03)0.52(0.04)
2000.64(0.03)0.66(0.08)
High500.61(0.03)0.58(0.03)
2000.67(0.04)0.68(0.03)
NormalbLow500.67(0.04)0.57(0.06)
2000.81(0.02)0.75(0.04)
High500.86(0.06)0.69(0.07)
2000.9(0.03)0.81(0.04)
BinarybLow500.57(0.06)0.51(0.09)
2000.63(0.05)0.6(0.1)
High500.68(0.06)0.6(0.07)
2000.78(0.02)0.76(0.03)
NormalcLow500.66(0.04)0.57(0.07)
2000.79(0.03)0.75(0.03)
High500.91(0.03)0.8(0.05)
2000.96(0.01)0.88(0.04)
BinarycLow500.59(0.05)0.5(0.04)
2000.73(0.04)0.67(0.04)
High500.94(0.03)0.81(0.04)
2000.97(0.01)0.92(0.01)

  1. 10% relevant (Q=60).

  2. 25% relevant (Q=150)

  3. 50% relevant (Q=300).

Table 10

(Simulation Study: Prediction Accuracy) Summary measures for AUCs obtained from true case-control status and BVS CLR/lasso CLR predicted case probabilities in independent test sets for K=20 and K=100 features (β1,…,βQ between 0.3 and 0.7 in magnitude).

Scenario 1Scenario 2
20 features100 features
Median(IQR)Median(IQR)
Data typeCorrelation levelNumber of pairsBVS CLRLasso CLRBVS CLRLasso CLR
NormalaLow500.45(0.19)0.5(0.02)0.54(0.17)0.5(0.08)
2000.48(0.14)0.46(0.15)0.58(0.15)0.56(0.16)
High500.49(0.24)0.5(<0.01)0.45(0.1)0.48(0.12)
2000.49(0.1)0.46(0.14)0.52(0.1)0.48(0.09)
BinaryaLow500.5(0.11)0.5(<0.01)0.52(0.09)0.5(0.01)
2000.47(0.13)0.48(0.13)0.52(0.1)0.54(0.16)
High500.49(0.11)0.5(<0.01)0.52(0.1)0.5(0.01)
2000.47(0.12)0.48(0.11)0.5(0.1)0.53(0.16)
NormalbLow500.79(0.17)0.8(0.25)0.65(0.19)0.5(0.13)
2000.85(0.05)0.86(0.05)0.68(0.19)0.67(0.22)
High500.66(0.21)0.71(0.3)0.87(0.04)0.87(0.04)
2000.71(0.1)0.74(0.09)0.88(0.03)0.89(0.02)
BinarybLow500.51(0.14)0.5(<0.01)0.56(0.12)0.5(0.08)
2000.52(0.14)0.52(0.17)0.58(0.09)0.58(0.15)
High500.59(0.14)0.5(0.06)0.54(0.1)0.5(0.09)
2000.64(0.11)0.65(0.11)0.57(0.09)0.53(0.1)
NormalcLow500.65(0.23)0.65(0.28)0.71(0.17)0.5(0.14)
2000.63(0.18)0.64(0.18)0.77(0.18)0.75(0.23)
High500.62(0.18)0.5(0.15)0.9(0.02)0.89(0.02)
2000.64(0.1)0.6(0.1)0.92(0.01)0.91(0.01)
BinarycLow500.62(0.15)0.5(0.18)0.58(0.12)0.5(0.13)
2000.66(0.09)0.66(0.1)0.63(0.11)0.6(0.16)
High500.6(0.13)0.5(0.14)0.5(0.08)0.46(0.09)
2000.65(0.1)0.66(0.09)0.55(0.07)0.51(0.08)

  1. 10% relevant (Q=2 for K=20; Q=10 for K=100).

  2. 25% relevant (Q=5 for K=20; Q=25 for K=100).

  3. 50% relevant (Q=10 for K=20; Q=50 for K=100).

Table 11

(Simulation Study: MSE for relevant predictors) Summary measures for BVS CLR/lasso CLR MSE estimates for sets of relevant predictors computed across predictors (β1,…,βQ at least 1 in magnitude).

Scenario 1Scenario 2Scenario 3
20 features100 features600 features
MedianMedianMedian
Data typeCorrelation levelRelevantNumber of pairsBVS CLRLasso CLRBVS CLRLasso CLRBVS CLRLasso CLR
NormalaLowYes500.290.220.461.021.151.30
2000.060.070.470.310.811.26
HighYes500.190.530.781.111.171.28
2000.050.100.750.961.031.24
BinaryaLowYes500.560.841.011.401.922.16
2000.090.170.290.401.341.83
HighYes500.650.990.841.112.082.13
2000.100.180.320.511.491.85
NormalbLowYes500.270.530.831.281.261.34
2000.070.120.270.701.071.30
HighYes500.390.601.091.251.171.31
2000.200.261.101.1811.28
BinarybLowYes500.730.761.681.972.132.25
2000.130.180.640.781.802.12
HighYes501.111.221.832.082.072.23
2000.140.260.630.781.702.05
NormalcLowYes500.470.760.921.241.251.31
2000.070.260.530.921.121.30
HighYes500.561.031.071.291.181.31
2000.400.880.961.191.071.30
BinarycLowYes500.970.901.942.172.202.24
2000.260.211.201.241.902.20
HighYes500.740.801.711.902.042.18
2000.260.191.041.131.882.10

  1. 10% relevant (Q=2 for K=20; Q=10 for K=100; Q=60 for K=600).

  2. 25% relevant (Q=5 for K=20; Q=25 for K=100; Q=150 for K=600).

  3. 50% relevant (Q=10 for K=20; Q=50 for K=100; Q=300 for K=600).

Table 12

(Simulation Study: MSE for relevant predictors) Summary measures for BVS CLR/lasso CLR MSE estimates for sets of relevant predictors computed across predictors (β1,…,βQ between 0.3 and 0.7 in magnitude).

Scenario 1Scenario 2
20 features100 features
MedianMedian
Data typeCorrelation levelRelevantNumber of pairsBVS CLRLasso CLRBVS CLRLasso CLR
NormalaLowYes500.160.140.180.19
2000.030.030.720.06
HighYes500.180.270.180.17
2000.020.060.210.14
BinaryaLowYes500.150.180.150.20
2000.100.110.230.12
HighYes500.190.220.170.20
2000.110.120.220.12
NormalbLowYes500.250.140.150.25
2000.040.030.460.09
HighYes500.140.160.180.20
2000.100.120.300.18
BinarybLowYes500.160.150.160.21
2000.100.110.240.13
HighYes500.220.270.170.21
2000.140.130.220.13
NormalcLowYes500.200.130.200.29
2000.070.030.240.15
HighYes500.190.200.180.24
2000.130.070.160.22
BinarycLowYes500.190.230.200.25
2000.130.100.200.14
HighYes500.250.310.210.24
2000.140.110.200.19

  1. 10% relevant (Q=2 for K=20; Q=10 for K=100).

  2. 25% relevant (Q=5 for K=20; Q=25 for K=100);

  3. 50% relevant (Q=10 for K=20; Q=50 for K=100).

Table 13

(Simulation Study: MSE for non-relevant predictors) Summary measures for BVS CLR/lasso CLR MSE estimates for sets of non-relevant predictors computed across predictors (β1,…,βQ at least 1 in magnitude).

Scenario 1Scenario 2Scenario 3
20 features100 features600 features
MedianMedianMedian
Data typeCorrelation levelRelevantNumber of pairsBVS CLRLasso CLRBVS CLRLasso CLRBVS CLRLasso CLR
NormalaLowNo500.030.010.02<0.010.02<0.01
200<0.01<0.010.01<0.010.04<0.01
HighNo500.030.010.03<0.010.01<0.01
200<0.010.010.03<0.010.01<0.01
BinaryaLowNo500.040.020.040.010.01<0.01
2000.010.010.030.010.03<0.01
HighNo500.040.030.040.010.01<0.01
2000.010.010.030.010.03<0.01
NormalbLowNo500.020.010.02<0.010.01<0.01
2000.010.010.01<0.010.04<0.01
HighNo500.020.010.01<0.010.01<0.01
200<0.010.010.01<0.010.01<0.01
BinarybLowNo500.040.060.020.010.01<0.01
2000.010.020.020.020.03<0.01
HighNo500.030.050.040.010.01<0.01
2000.010.030.010.020.02<0.01
NormalcLowNo500.020.010.02<0.010.02<0.01
200<0.010.010.01<0.010.04<0.01
HighNo500.010.010.01<0.01<0.01<0.01
200<0.010.010.01<0.010.01<0.01
BinarycLowNo500.030.060.030.010.02<0.01
2000.010.040.010.020.03<0.01
HighNo500.040.070.040.010.01<0.01
2000.010.040.010.020.01<0.01

  1. 10% relevant (Q=2 for K=20; Q=10 for K=100; Q=60 for K=600).

  2. 25% relevant (Q=5 for K=20; Q=25 for K=100; Q=150 for K=600);

  3. 50% relevant (Q=10 for K=20; Q=50 for K=100; Q=300 for K=600).

Table 14

(Simulation Study: MSE for non-relevant predictors) Summary measures for BVS CLR/lasso CLR MSE estimates for sets of relevant predictors computed across predictors (β1,…,βQ between 0.3 and 0.7 in magnitude).

Scenario 1Scenario 2
20 features100 features
MedianMedian
Data typeCorrealtion levelRelevantNumber of pairsBVS CLRLasso CLRBVS CLRLasso CLR
NormalaLowYes500.03<0.010.04<0.01
200<0.01<0.010.07<0.01
HighYes500.04<0.010.05<0.01
200<0.01<0.010.04<0.01
BinaryaLowYes500.050.010.05<0.01
2000.010.010.04<0.01
HighYes500.050.010.05<0.01
2000.01<0.010.04<0.01
NormalbLowYes500.040.010.03<0.01
200<0.010.010.04<0.01
HighYes500.040.010.02<0.01
200<0.010.010.04<0.01
BinarybLowYes500.050.010.05<0.01
2000.010.010.050.01
HighYes500.050.010.050.01
2000.010.010.040.01
NormalcLowYes500.040.020.03<0.01
200<0.010.010.030.01
HighYes500.030.010.01<0.01
200<0.010.010.020.01
BinarycLowYes500.040.020.040.01
2000.010.020.040.01
HighYes500.040.020.040.01
2000.010.020.040.01

  1. 10% relevant (Q=2 for K=20; Q=10 for K=100).

  2. 25% relevant (Q=5 for K=20; Q=25 for K=100).

  3. 50% relevant (Q=10 for K=20; Q=50 for K=100).

B Trace Plots

Figure 4: Log posterior probability plots for simulation settings with K$$K$$ normally distributed features of which Q$$Q$$ are relevant and weakly correlated (I=50$$I=50$$ pairs).
Figure 4:

Log posterior probability plots for simulation settings with K normally distributed features of which Q are relevant and weakly correlated (I=50 pairs).

Figure 5: Number of selected features across post burn-in iterations for simulation settings with K$$K$$ normally distributed features of which Q$$Q$$ are relevant and weakly correlated (I=50$$I=50$$ pairs).
Figure 5:

Number of selected features across post burn-in iterations for simulation settings with K normally distributed features of which Q are relevant and weakly correlated (I=50 pairs).

Figure 6: Inclusion status (0 = no, 1 = yes) across burn-in iterations for features with top 3 inclusion probabilities for simulation settings with K = 20 normally distributed features of which Q = 5 are relevant and weakly correlated (I=50$$I=50$$ pairs).
Figure 6:

Inclusion status (0 = no, 1 = yes) across burn-in iterations for features with top 3 inclusion probabilities for simulation settings with K = 20 normally distributed features of which Q = 5 are relevant and weakly correlated (I=50 pairs).

References

1. Balasubramanian R, Houseman EA, Coull BA, Lev M, Schwamm LH, Betensky RA. Variable importance in matched case-control studies in settings of high dimensional data. J Roy Stat Soc Ser C 2014;63:639–655.10.1111/rssc.12056Search in Google Scholar

2. Anglim PP, Galler JS, Koss MN, Hagen JA, Turla S, Campan M, et al. Identification of a panel of sensitive and specific DNA methylation markers for squamous cell lung cancer. Mol Cancer 2008;7:62.10.1186/1476-4598-7-62Search in Google Scholar PubMed PubMed Central

3. Westman E, Simmons A, Zhang Y, Muehlboeck JS, Tunnard C, Liu Y, et al. Multivariate analysis of MRI data for Alzheimer’s disease, mild cognitive impairment and healthy controls. Neuroimage 2011;54:1178–1187.10.1016/j.neuroimage.2010.08.044Search in Google Scholar PubMed

4. Breslow N, Day N. Statistical methods in cancer research (vol. 1): the analysisof case-control studies. Lyon: IARC, 1980.Search in Google Scholar

5. Tan Q, Thomassen M, Kruse TA. Feature selection for predicting tumor metastases in microarray experiments using paired design. Cancer Inf 2007;3:213–218.10.1177/117693510700300025Search in Google Scholar

6. Adewale AJ, Dinu I, Yasui Y. Boosting for correlated binary classification. J Comput Graph Stat 2010;19:140–153.10.1198/jcgs.2009.07118Search in Google Scholar

7. Friedman J. Greedy function approximation: a gradient boosting machine. Ann Stat 2001;29:1189–1232.10.1214/aos/1013203451Search in Google Scholar

8. Freund Y, Schapire R. A short introduction to boosting. J Jpn Soc Artif Intell 1999;14:771–780.Search in Google Scholar

9. Frank I, Friedman J. A statistical view of some chemometrics regression tools. Technometrics 1993;35:109–135.10.1080/00401706.1993.10485033Search in Google Scholar

10. Breiman L. Random forests. Mach Learn 2001;45:5–32.10.1023/A:1010933404324Search in Google Scholar

11. Qian J, Payabvash S, Kemmling A, Lev M, Schwamm L, Betensky RA. Variable selection and prediction using a nested, matched case-control study: application to hospital acquired pneumonia in stroke patients. Biometrics, 2013; 70:153-163. doi:10.1111/biom.12113.Search in Google Scholar PubMed PubMed Central

12. Tibshirani R. Regression shrinkage and selection via the lasso. J Roy Stat Soc Ser B Methodol 1996;58:267–288.10.1111/j.2517-6161.1996.tb02080.xSearch in Google Scholar

13. Zou H, Hastie T. Regularization and variable selection via the elastic net. J Roy Stat Soc Ser B 2005;67:301–320.10.1111/j.1467-9868.2005.00503.xSearch in Google Scholar

14. Mitchell T, Beauchamp J. Bayesian variable selection in linear regression. J Am Stat Assoc 1988;83:1023–1032.10.1080/01621459.1988.10478694Search in Google Scholar

15. George E, McCulloch R. Variable selection via Gibbs sampling. J Am Stat Assoc 1993;88:881–889.10.1080/01621459.1993.10476353Search in Google Scholar

16. Ishwaran H, Rao JS. Spike and slab variable selection: frequentist and Bayesian strategies. Ann Stat 2005;33:730–773.10.1214/009053604000001147Search in Google Scholar

17. George E, McCulloch R. Approaches for Bayesian variable selection. Stat Sin 1997;7:339–373.Search in Google Scholar

18. Barbieri M, Berger J. Optimal predictive model selection. Ann Stat 2004;32:870–897.10.1214/009053604000000238Search in Google Scholar

19. Ishwaran H, Rao S. Clustering gene expression profile data by selective shrinkage. Stat Probab Lett 2008;78:1490–1497.10.1016/j.spl.2008.01.003Search in Google Scholar

20. Rockova V, George E. The spike-and-slab LASSO. Submitted manuscript, 2015:1–39.Search in Google Scholar

21. Rockova V. Bayesian estimation of sparse signals with a continuous spike-and-slab prior. Ann Stat 2015:1–44. (under revision).Search in Google Scholar

22. Pang X, Gill J. Spike and slab prior distributions for simultaneous Bayesian hypothesis testing, model selection, and prediction, of nonlinear outcomes. Washington University in St. Louis 2009.Search in Google Scholar

23. Lee KE, Sha N, Dougherty ER, Vannucci M, Mallick BK. Gene selection: a Bayesian variable selection approach. Bioinformatics 2003;19:90–97.10.1093/bioinformatics/19.1.90Search in Google Scholar PubMed

24. Sha N, Vannucci M, Tadesse MG, Brown PJ, Dragoni I, Davies N, et al. Bayesian variable selection in multinomial probit models to identify molecular signatures of disease stage. Biometrics 2004;60:812–819.10.1111/j.0006-341X.2004.00233.xSearch in Google Scholar PubMed

25. Zhou X, Liu KY, Wong ST. Cancer classification and prediction using logistic regression with Bayesian gene selection. J Biomed Inf 2004;37:249–259.10.1016/j.jbi.2004.07.009Search in Google Scholar PubMed

26. Celeux G, Anbari ME, Marin JM, Robert CP. Regularization in regression: comparing Bayesian and frequentist methods in a poorly informative situation. Bayesian Anal 2012;7:477–502.10.1214/12-BA716Search in Google Scholar

27. Johnson V. On numerical aspects of Bayesian model selection in high and ultrahigh-dimensional settings. Bayesian Anal 2004;1:1–17.10.1214/13-BA818Search in Google Scholar PubMed PubMed Central

28. Smith M, Fahrmeir L. Spatial Bayesian variable selection with application to functional magnetic resonance imaging. J Am Stat Assoc 2007;102:417–431.10.1198/016214506000001031Search in Google Scholar

29. Stingo F, Chen Y, Tadesse M, Vannucci M. Incorporating biological information into linear models: a Bayesian approach to the selection of pathways and genes. Ann Appl Stat 2011;5:1978–2002.10.1214/11-AOAS463Search in Google Scholar PubMed PubMed Central

30. Yuan M, Lin Y. Model selection and estimation in regression with grouped variables. J Roy Stat Soc Ser B Stat Methodol 2006;68:49–67.10.1111/j.1467-9868.2005.00532.xSearch in Google Scholar

31. Meier L, van de Geer S, Bühlmann P. The group lasso for logistic regression. J Roy Stat Soc Ser B Stat Methodol 2008;70:53–71.10.1111/j.1467-9868.2007.00627.xSearch in Google Scholar

32. Friedman J, Hastie T, Tibshirani R. A note on the group lasso and a sparse group lasso. Technical report. Department of Statistics, Stanford University, 2010.Search in Google Scholar

33. Fan J, Li R. Variable selection via nonconcave penalized likelihood and its Oracle properties. J Am Stat Assoc 2001;96:1348–1360.10.1198/016214501753382273Search in Google Scholar

34. Zou H. The adaptive lasso and its Oracle properties. J Am Stat Assoc 2006;101:1418–1429.10.1198/016214506000000735Search in Google Scholar

35. Gelman A, Roberts G, Gilks W. Efficient metropolis jumping rules. Bayesian Stat 1996;5:599–607.Search in Google Scholar

36. Roberts G, Gelman A, Gilks W. Weak convergence and optimal scaling of random walk metropolis algorithms. Ann Appl Probab 1997;7:110–120.Search in Google Scholar

37. Lamnisos D, Griffin JE, Steel MF. Adaptive Monte Carlo for Bayesian variable selection in regression models. J Comput Graph Stat 2013;22:729–748. doi:10.1080/10618600.2012.694756.Search in Google Scholar

38. Gelman A, Rubin D. Inference from iterative simulation using multiple sequences. Stat Sci 1992;7:457–511.10.1214/ss/1177011136Search in Google Scholar

39. Brooks S, Gelman A. General methods for monitoring convergence of iterative simulations. J Comput Graph Stat 1998;7:434–455.Search in Google Scholar

40. Javanmard A, Montanari A. Confidence intervals and hypothesis testing for high-dimensional regression. J Mach Learn Res 2014;15:2869–2909.Search in Google Scholar

41. Reid S, Tibshirani R, Friedman J. A study of error variance estimation in lasso regression. Stat Sin 26 (2016), 35-67. doi:10.5705/ss.2014.042.Search in Google Scholar

42. Core Team R, Language R: A. and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria, 2014;Available at http://www.R-project.org/.Search in Google Scholar

43. Kemmling A, Lev M, Payabvash S, Betensky R, Qian J, Masrur S, et al. Hospital acquired pneumonia is linked to right peri-insular stroke. PLoS ONE 2013;8:e71141. doi:10.1371/journal.pone.0071141.Search in Google Scholar PubMed PubMed Central

44. Desikan R, Ségonne F, Fischl B, Quinn B, Dickerson B, Blacker D, et al. An automated labeling system for subdividing the human cerebral cortex on MRI scans into gyral based regions of interest. Neuroimage 2006;31:968–980.10.1016/j.neuroimage.2006.01.021Search in Google Scholar PubMed

45. Hua K, Zhang J, Wakana S, Jiang H, Li X, Reich D, et al. Tract probability maps in stereotaxic spaces: analyses of white matter anatomy and tract-specific quantification. Neuroimage 2008;39:336–347.10.1016/j.neuroimage.2007.07.053Search in Google Scholar PubMed PubMed Central

46. Cechetto D, Chen S. Subcortical sites mediating sympathetic responses from insular cortex in rats. Am J Physiol 1990;258:R245–255.10.1152/ajpregu.1990.258.1.R245Search in Google Scholar PubMed

47. Sander D, Klingelhöfer J. Changes of circadian blood pressure patterns and cardiovascular parameters indicate lateralization of sympathetic activation following hemispheric brain infarction. J Neurol 1995;242:313–318.10.1007/BF00878874Search in Google Scholar PubMed

48. Meyer S, Strittmatter M, Fischer C, Georg T, Schmitz B. Lateralization in autononic dysfunction in ischemic stroke involving the insular cortex. Neuroreport 2004;15:357–361.10.1097/00001756-200402090-00029Search in Google Scholar PubMed

49. Colivicchi F, Bassi A, Santini M, Caltagirone C. Cardiac autonomic derangement and arrhythmias in right-sided stroke with insular involvement. Stroke J Cereb Circ 2004;35:2094–2098.10.1161/01.STR.0000138452.81003.4cSearch in Google Scholar PubMed

50. Muller P, Parmigiani G, FDR Rice K. Bayesian multiple comparisons rules. Technical report. Johns Hopkins University, Dept. of Biostatistics Working Papers. Working Paper 115. 2006;Available at http://biostats.bepress.com/jhubiostat/paper115.Search in Google Scholar

51. Cassese A, Guindani M, Tadesse MG, Falciani F, Vannucci M. A hierarchical Bayesian model for inference of copy number variants and their association to gene expression. Ann Appl Stat 2014;8:148–175. doi:10.1214/13-AOAS705.Search in Google Scholar PubMed PubMed Central

52. Lewin A, Saadi H, Peters JE, Moreno-Moral A, Lee JC, Smith KG, et al. MT-HESS: an efficient Bayesian approach for simultaneous association detection in OMICS datasets, with application to eQTL mapping in multiple tissues. Bioinformatics 2016;32:523–532. doi:10.1093/bioinformatics/btv568.Search in Google Scholar PubMed PubMed Central

53. Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J Roy Stat Soc Ser B 1995;57:289–300.10.1111/j.2517-6161.1995.tb02031.xSearch in Google Scholar

54. Guo Y , Balasubramanian R. Comparative evaluation of classifiers in the presence of statistical interactions between features in high dimensional data settings. Int J Biostat 2012;8. Article 17. doi:10.1515/1557-4679.1373.Search in Google Scholar PubMed

55. Chipman H. Bayesian variable selection with related predictors. Canadian J Stat 1996;24:17–36.10.2307/3315687Search in Google Scholar

56. Mukherjee B, Liu I, Sinha S. Analysis of matched case-control data with multiple ordered disease states: possible choices and comparisons. Stat Med 2007;26:3240–3257.10.1002/sim.2790Search in Google Scholar PubMed

57. Carroll R, Ruppert D, Stefanski L. Monographs on statistics and applied probability. Measurement error in nonlinear models, 1st ed. Vol. 63. Boca Raton: Chapman & Hall/ CRC, 1995.10.1007/978-1-4899-4477-1Search in Google Scholar

Published Online: 2017-1-31

© 2017 Walter de Gruyter GmbH, Berlin/Boston

Articles in the same Issue

  1. Commentary
  2. Big Data, Small Sample
  3. Research Articles
  4. Parameter Estimation of a Two-Colored Urn Model Class
  5. Combinatorial Mixtures of Multiparameter Distributions: An Application to Bivariate Data
  6. On the Conditional Power in Survival Time Analysis Considering Cure Fractions
  7. Comparing Four Methods for Estimating Tree-Based Treatment Regimes
  8. On Stratified Adjusted Tests by Binomial Trials
  9. Improvement Screening for Ultra-High Dimensional Data with Censored Survival Outcomes and Varying Coefficients
  10. Bayesian Variable Selection Methods for Matched Case-Control Studies
  11. Testing Equality of Treatments under an Incomplete Block Crossover Design with Ordinal Responses
  12. Empirical Likelihood in Nonignorable Covariate-Missing Data Problems
  13. A Quantitative Concordance Measure for Comparing and Combining Treatment Selection Markers
  14. Median Analysis of Repeated Measures Associated with Recurrent Events in Presence of Terminal Event
  15. A Theorem at the Core of Colliding Bias
  16. Group Tests for High-dimensional Failure Time Data with the Additive Hazards Models
  17. Characterizing Highly Benefited Patients in Randomized Clinical Trials
https://doi.org/10.1515/ijb-2016-0043
Keywords for this article
Bayesian analysis; conditional logistic regression; matched case-control studies; variable selection methods

Articles in the same Issue

  1. Commentary
  2. Big Data, Small Sample
  3. Research Articles
  4. Parameter Estimation of a Two-Colored Urn Model Class
  5. Combinatorial Mixtures of Multiparameter Distributions: An Application to Bivariate Data
  6. On the Conditional Power in Survival Time Analysis Considering Cure Fractions
  7. Comparing Four Methods for Estimating Tree-Based Treatment Regimes
  8. On Stratified Adjusted Tests by Binomial Trials
  9. Improvement Screening for Ultra-High Dimensional Data with Censored Survival Outcomes and Varying Coefficients
  10. Bayesian Variable Selection Methods for Matched Case-Control Studies
  11. Testing Equality of Treatments under an Incomplete Block Crossover Design with Ordinal Responses
  12. Empirical Likelihood in Nonignorable Covariate-Missing Data Problems
  13. A Quantitative Concordance Measure for Comparing and Combining Treatment Selection Markers
  14. Median Analysis of Repeated Measures Associated with Recurrent Events in Presence of Terminal Event
  15. A Theorem at the Core of Colliding Bias
  16. Group Tests for High-dimensional Failure Time Data with the Additive Hazards Models
  17. Characterizing Highly Benefited Patients in Randomized Clinical Trials
Sign up now to receive a 20% welcome discount
Institutional Access
How does access work?
Have an idea on how to improve our website?
Please write us.
© 2025 De Gruyter Brill
Downloaded on 24.9.2025 from https://www.degruyterbrill.com/document/doi/10.1515/ijb-2016-0043/html
Scroll to top button