Resampling-based multiple testing methods with covariate adjustment: application to investigation of antiretroviral drug susceptibility

Summary .   Identifying genetic mutations that cause clinical resistance to antiretroviral drugs requires adjustment for potential confounders, such as the number of active drugs in a HIV-infected patient's regimen other than the one of interest. Motivated by this problem, we investigated resampling-based methods to test equal mean response across multiple groups defined by HIV genotype, after adjustment for covariates. We consider construction of test statistics and their null distributions under two types of model: parametric and semiparametric. The covariate function is explicitly specified in the parametric but not in the semiparametric approach. The parametric approach is more precise when models are correctly specified, but suffer from bias when they are not; the semiparametric approach is more robust to model misspecification, but may be less efficient. To help preserve type I error while also improving power in both approaches, we propose resampling approaches based on matching of observations with similar covariate values. Matching reduces the impact of model misspecification as well as imprecision in estimation. These methods are evaluated via simulation studies and applied to a data set that combines results from a variety of clinical studies of salvage regimens. Our focus is on relating HIV genotype to viral susceptibility to abacavir after adjustment for the number of active antiretroviral drugs (excluding abacavir) in the patient's regimen.     

Information Sets and Excess Zeros in Random Effects Modeling of Longitudinal Data

Marginal regression via generalized estimating equations is widely used in biostatistics to model longitudinal data from subjects whose outcomes and covariates are observed at several time points. In this paper we consider two issues that have been raised in the literature concerning the marginal regression approach. The first is that even though the past history may be predictive of outcome, the marginal approach does not use this history. Although marginal regression has the flexibility of allowing between-subject variations in the observation times, it may lose substantial prediction power in comparison with the transitional modeling approach that relates the responses to the covariate and outcome histories. We address this issue by using the concept of “information sets” for prediction to generalize the “partly conditional mean” approach of Pepe and Couper (J. Am. Stat. Assoc. 92:991–998, 1997). This modeling approach strikes a balance between the flexibility of the marginal approach and the predictive power of transitional modeling. Another issue is the problem of excess zeros in the outcomes over what the underlying model for marginal regression implies. We show how our predictive modeling approach based on information sets can be readily modified to handle the excess zeros in the longitudinal time series. By synthesizing the marginal, transitional, and mixed effects modeling approaches in a predictive framework, we also discuss how their respective advantages can be retained while their limitations can be circumvented for modeling longitudinal data.     

A Comparison Between Several Goodness of Fit Tests for Cox's Model

A comparison is made between two approaches to testing goodness of fit of Cox's regression model for survival data. The first approach is based on the inclusion of time dependent covariates, whereas the second one is based on the autocovariance of successive contributions to the derivative of the loglikelihood. It appears that the second test is most appropriate for testing in situations where the structure of the departure from proportional hazards is not known a priori. An approximate expression for the relative efficiency of the two test procedures is presented.     

Use of weighted reference panels based on empirical estimates of ancestry for capturing untyped variation

Many association methods use a subset of genotyped single nucleotide polymorphisms (SNPs) to capture or infer genotypes at other untyped SNPs. We and others previously showed that tag SNPs selected to capture common variation using data from The International HapMap Consortium (Nature 437:1299–1320, 2005), The International HapMap Consortium (Nature 449:851–861, 2007) could also capture variation in populations of similar ancestry to HapMap reference populations (de Bakker et al. in Nat Genet 38:1298–1303, 2006; González-Neira et al. in Genome Res 16:323–330, 2006; Montpetit et al. in PLoS Genet 2:282–290, 2006; Mueller et al. in Am J Hum Genet 76:387–398, 2005). To capture variation in admixed populations or populations less similar to HapMap panels, a “cosmopolitan approach,” in which all samples from HapMap are used as a single reference panel, was proposed. Here we refine this suggestion and show that use of a “weighted reference panel,” constructed based on empirical estimates of ancestry in the target population (relative to available reference panels), is more efficient than the cosmopolitan approach. Weighted reference panels capture, on average, only slightly fewer common variants (minor allele frequency > 5%) than the cosmopolitan approach (mean r ² = 0.977 vs. 0.989, 94.5% variation captured vs. 96.8% at r ² > 0.8), across the five populations of the Multiethnic Cohort, but entail approximately 25% fewer tag SNPs per panel (average 538 vs. 718). These results extend a recent study in two Indian populations (Pemberton et al. in Ann Hum Genet 72:535–546, 2008). Weighted reference panels are potentially useful for both the selection of tag SNPs in diverse populations and perhaps in the design of reference panels for imputation of untyped genotypes in genome-wide association studies in admixed populations. Electronic supplementary material  The online version of this article (doi:) contains supplementary material, which is available to authorized users.     

A unique genome-wide association analysis in extended Utah high-risk pedigrees identifies a novel melanoma risk variant on chromosome arm 10q

Only two genome-wide association (GWA) screens have been published for melanoma (Nat Genet 47:920–925, 2009; Nat Genet 40:838–840, 2008). Using a unique approach, we performed a genome-wide association study in 156 related melanoma cases from 34 high-risk Utah pedigrees. Genome-wide association analysis was performed on nearly 500,000 markers; we compared cases to 2,150 genotypically matched samples from Illumina’s iControls database. We performed genome-wide association with EMMAX software, which is designed to account for population structure, including relatedness between cases. Three SNPs exceeded a genome-wide significance threshold of p < 5 × 10⁻⁸ on chromosome arm 10q25.1 (rs17119434, rs17119461, and rs17119490), where the most extreme p value was 7.21 × 10⁻¹². This study represents a new and unique approach to predisposition gene identification; and it is the first genome-wide association study performed in related cases in high-risk pedigrees. Our approach illustrates an example of using high-risk pedigrees for the identification of new melanoma predisposition variants.     

A Propensity Score-Enhanced Sequential Analytic Method for Comparative Drug Safety Surveillance

We introduce a new sequential monitoring approach to facilitate the use of observational electronic healthcare utilization databases in comparative drug safety surveillance studies comparing the safety between two approved medical products. The new approach enhances the confounder adjustment capabilities of the conditional sequential sampling procedure (CSSP), an existing group sequential method for sequentially monitoring excess risks of adverse events following the introduction of a new medical product. It applies to a prospective cohort setting where information for both treatment and comparison groups accumulates concurrently over time. CSSP adjusts for covariates through stratification and thus it may have limited capacity to control for confounding as it can only accommodate a few categorical covariates. To address this issue, we propose the propensity score (PS)-stratified CSSP, in which we construct strata based on selected percentiles of the estimated PSs. The PS is defined as the conditional probability of being treated given measured baseline covariates and is commonly used in epidemiological studies to adjust for confounding bias. The PS-stratified CSSP approach integrates this more flexible confounding adjustment, PS-stratification, with the sequential analytic approach, CSSP, thus inheriting CSSP’s attractive features: (i) it accommodates varying amounts of person follow-up time, (ii) it uses exact conditional inference, which can be important when studying rare safety outcomes, and (iii) it allows for a large number of interim tests. Further, it overcomes CSSP’s difficulty with adjusting for multiple categorical and continuous confounders.     

Parameter estimation and goodness-of-fit in log binomial regression

An estimate of the risk, adjusted for confounders, can be obtained from a fitted logistic regression model, but it substantially over-estimates when the outcome is not rare. The log binomial model, binomial errors and log link, is increasingly being used for this purpose. However this model's performance, goodness of fit tests and case-wise diagnostics have not been studied. Extensive simulations are used to compare the performance of the log binomial, a logistic regression based method proposed by Schouten et al. (1993) and a Poisson regression approach proposed by Zou (2004) and Carter, Lipsitz, and Tilley (2005). Log binomial regression resulted in "failure" rates (non-convergence, out-of-bounds predicted probabilities) as high as 59%. Estimates by the method of Schouten et al. (1993) produced fitted log binomial probabilities greater than unity in up to 19% of samples to which a log binomial model had been successfully fit and in up to 78% of samples when the log binomial model fit failed. Similar percentages were observed for the Poisson regression approach. Coefficient and standard error estimates from the three models were similar. Rejection rates for goodness of fit tests for log binomial fit were around 5%. Power of goodness of fit tests was modest when an incorrect logistic regression model was fit. Examples demonstrate the use of the methods. Uncritical use of the log binomial regression model is not recommended.     

Covariance adjustment of survival curves based on Cox's proportional hazards regression model

Cox^'s proportional hazards regression model is a useful statisticaltool for the analysis of ‘survival data’ from longitudinalstudies. This multivariate method compares the ‘survivalexperience’ between two or more exposure groups whileallowing for simultaneous adjustment of confounding due to oneor more covariates. In addition to the summary regression statistics,further insight on the exposure–response relationshipcan be gained by visually examining the covariates–adjustedsurvival curves in the respective comparison groups. Covariates–adjustedsurvival curves are usually computed by the ‘average covariatemethod’. This method is, however, subject to potentialdrawbacks. A method that avoids these drawbacks is to estimateadjusted survival curves by the corrected group prognostic curvesapproach. We have written a computer program to construct survivalcurves by the latter method. The program is coded in the InteractiveMatrix Language of SAS.     

Testing the productive-space hypothesis: rational and power

Understanding and explaining the causes of variation in food-chain length is a fundamental challenge for community ecology. The productive-space hypothesis, which suggests food-chain length is determined by the combination of local resource availability and ecosystem size, is central to this challenge. Two different approaches currently exist for testing the productive-space hypothesis: (1) the dual gradient approach that tests for significant relationships between food-chain length and separate gradients of ecosystem size (e.g., lake volume) and per-unit-size resource availability (e.g., g C m⁻¹ year⁻²), and (2) the single gradient approach that tests for a significant relationship between food-chain length and the productive space (product of ecosystem size and per-unit-size resource availability). Here I evaluate the efficacy of the two approaches for testing the productive-space hypothesis. Using simulated data sets, I estimate the Type 1 and Type 2 error rates for single and dual gradient models in recovering a known relationship between food-chain length and ecosystem size, resource availability, or the combination of ecosystem size and resource ability, as specified by the productive-space hypothesis. The single gradient model provided high power (low Type 2 error rates) but had a very high Type 1 error rate, often erroneously supporting the productive-space hypothesis. The dual gradient model had a very low Type 1 error rate but suffered from low power to detect an effect of per-unit-size resource availability because the range of variation in resource availability is limited. Finally, I performed a retrospective power analysis for the Post et al. (Nature 405:1047–1049, 2000) data set, which tested and rejected the productive-space hypothesis using the dual gradient approach. I found that Post et al. (Nature 405:1047–1049, 2000) had sufficient power to reject the productive-space hypothesis in north temperate lakes; however, the productive-space hypothesis must be tested in other ecosystems before its generality can be fully addressed.     

Estimating differences in restricted mean lifetime using observational data subject to dependent censoring

In epidemiologic studies of time to an event, mean lifetime is often of direct interest. We propose methods to estimate group- (e.g., treatment-) specific differences in restricted mean lifetime for studies where treatment is not randomized and lifetimes are subject to both dependent and independent censoring. The proposed methods may be viewed as a hybrid of two general approaches to accounting for confounders. Specifically, treatment-specific proportional hazards models are employed to account for baseline covariates, while inverse probability of censoring weighting is used to accommodate time-dependent predictors of censoring. The average causal effect is then obtained by averaging over differences in fitted values based on the proportional hazards models. Large-sample properties of the proposed estimators are derived and simulation studies are conducted to assess their finite-sample applicability. We apply the proposed methods to liver wait list mortality data from the Scientific Registry of Transplant Recipients.     

Perceptual learning with perceptions

In this work we present an approach to understand neuronal mechanisms underlying perceptual learning. Experimental results achieved with stimulus patterns of coherently moving dots are considered to build a simple neuronal model. The design of the model is made transparent and underlying behavioral assumptions made explicit. The key aspect of the suggested neuronal model is the learning algorithm used: We evaluated an implementation of Hebbian learning and are thus able to provide a straight-forward model capable to explain the neuronal dynamics underlying perceptual learning. Moreover, the simulation results suggest a very simple explanation for the aspect of “sub-threshold” learning (Watanabe et al. in Nature 413:844–884, 2001) as well as the relearning of motion discrimination after damage to primary visual cortex as recently reported (Huxlin et al. in J Neurosci 29:3981–3991, 2009) and at least indicate that perceptual learning might only occur when accompanied by conscious percepts.     

Background stratified Poisson regression analysis of cohort data

Background stratified Poisson regression is an approach that has been used in the analysis of data derived from a variety of epidemiologically important studies of radiation-exposed populations, including uranium miners, nuclear industry workers, and atomic bomb survivors. We describe a novel approach to fit Poisson regression models that adjust for a set of covariates through background stratification while directly estimating the radiation-disease association of primary interest. The approach makes use of an expression for the Poisson likelihood that treats the coefficients for stratum-specific indicator variables as ‘nuisance’ variables and avoids the need to explicitly estimate the coefficients for these stratum-specific parameters. Log-linear models, as well as other general relative rate models, are accommodated. This approach is illustrated using data from the Life Span Study of Japanese atomic bomb survivors and data from a study of underground uranium miners. The point estimate and confidence interval obtained from this ‘conditional’ regression approach are identical to the values obtained using unconditional Poisson regression with model terms for each background stratum. Moreover, it is shown that the proposed approach allows estimation of background stratified Poisson regression models of non-standard form, such as models that parameterize latency effects, as well as regression models in which the number of strata is large, thereby overcoming the limitations of previously available statistical software for fitting background stratified Poisson regression models.     

Dynamic formation of oriented patches in chondrocyte cell cultures

Growth factors have a significant impact not only on the growth dynamics but also on the phenotype of chondrocytes (Barbero et al. in J. Cell. Phys. 204:830–838, 2005). In particular, as chondrocytes approach confluence, the cells tend to align and form coherent patches. Starting from a mathematical model for fibroblast populations at equilibrium (Mogilner et al. in Physica D 89:346–367, 1996), a dynamic continuum model with logistic growth is developed. Both linear stability analysis and numerical solutions of the time-dependent nonlinear integro-partial differential equation are used to identify the key parameters that lead to pattern formation in the model. The numerical results are compared quantitatively to experimental data by extracting statistical information on orientation, density and patch size through Gabor filters.     

An evaluation of prior influence on the predictive ability of Bayesian model averaging

Model averaging is gaining popularity among ecologists for making inference and predictions. Methods for combining models include Bayesian model averaging (BMA) and Akaike’s Information Criterion (AIC) model averaging. BMA can be implemented with different prior model weights, including the Kullback–Leibler prior associated with AIC model averaging, but it is unclear how the prior model weight affects model results in a predictive context. Here, we implemented BMA using the Bayesian Information Criterion (BIC) approximation to Bayes factors for building predictive models of bird abundance and occurrence in the Chihuahuan Desert of New Mexico. We examined how model predictive ability differed across four prior model weights, and how averaged coefficient estimates, standard errors and coefficients’ posterior probabilities varied for 16 bird species. We also compared the predictive ability of BMA models to a best single-model approach. Overall, Occam’s prior of parsimony provided the best predictive models. In general, the Kullback–Leibler prior, however, favored complex models of lower predictive ability. BMA performed better than a best single-model approach independently of the prior model weight for 6 out of 16 species. For 6 other species, the choice of the prior model weight affected whether BMA was better than the best single-model approach. Our results demonstrate that parsimonious priors may be favorable over priors that favor complexity for making predictions. The approach we present has direct applications in ecology for better predicting patterns of species’ abundance and occurrence.     

Three algorithms and SAS macros for estimating power and sample size for logistic models with one or more independent variables of interest in the presence of covariates

Background

Commonly when designing studies, researchers propose to measure several independent variables in a regression model, a subset of which are identified as the main variables of interest while the rest are retained in a model as covariates or confounders. Power for linear regression in this setting can be calculated using SAS PROC POWER. There exists a void in estimating power for the logistic regression models in the same setting.

Methods

Currently, an approach that calculates power for only one variable of interest in the presence of other covariates for logistic regression is in common use and works well for this special case. In this paper we propose three related algorithms along with corresponding SAS macros that extend power estimation for one or more primary variables of interest in the presence of some confounders.

Results

The three proposed empirical algorithms employ likelihood ratio test to provide a user with either a power estimate for a given sample size, a quick sample size estimate for a given power, and an approximate power curve for a range of sample sizes. A user can specify odds ratios for a combination of binary, uniform and standard normal independent variables of interest, and or remaining covariates/confounders in the model, along with a correlation between variables.

Conclusions

These user friendly algorithms and macro tools are a promising solution that can fill the void for estimation of power for logistic regression when multiple independent variables are of interest, in the presence of additional covariates in the model.

     

A powerful and flexible multilocus association test for quantitative traits

Association mapping of complex traits typically employs tagSNP genotype data to identify a trait locus within a region of interest. However, considerable debate exists regarding the most powerful strategy for utilizing such tagSNP data for inference. A popular approach tests each tagSNP within the region individually, but such tests could lose power as a result of incomplete linkage disequilibrium between the genotyped tagSNP and the trait locus. Alternatively, one can jointly test all tagSNPs simultaneously within the region (by using genotypes or haplotypes), but such multivariate tests have large degrees of freedom that can also compromise power. Here, we consider a semiparametric model for quantitative-trait mapping that uses genetic information from multiple tagSNPs simultaneously in analysis but produces a test statistic with reduced degrees of freedom compared to existing multivariate approaches. We fit this model by using a dimension-reducing technique called least-squares kernel machines, which we show is identical to analysis using a specific linear mixed model (which we can fit by using standard software packages like SAS and R). Using simulated SNP data based on real data from the International HapMap Project, we demonstrate that our approach often has superior performance for association mapping of quantitative traits compared to the popular approach of single-tagSNP testing. Our approach is also flexible, because it allows easy modeling of covariates and, if interest exists, high-dimensional interactions among tagSNPs and environmental predictors.     

Optimal biocatalyst loading in a fixed bed

The optimal distribution of biocatalyst in a fixed bed operating at steady state was determined to minimize the length of the bed for a fixed conversion of 95%. The distribution in terms of the biocatalyst loading on an inert support depends upon the axial distance from the bed entrance (continuous solution) as well as a set of dimensionless parameters that reflect the bed geometry, the bulk flow, reaction kinetics and diffusional characteristics. A mathematical model of the system with the following features was used to obtain the results: (1) convective mass transfer and dispersion in the bulk phase; (2) mass transfer from the bulk phase to the surface of the catalyst particle; and (3) simultaneous diffusion and chemical reaction in the catalyst particle with Michaelis–Menton kinetics and a reliable diffusion model (Zhao and DeLancey in Biotechnol Bioeng 64:434–441, 1999, 2000). The solution to the mathematical model was obtained with Mathematica utilizing the Runge Kutta 4–5 method. The dimensionless length resulting from the continuous solution was compared with the optimal length restricted to a uniform or constant cell loading across the entire bed. The maximum difference in the dimensionless length between the continuous and uniform solutions was determined to be 6.5%. The model was applied to published conversion data for the continuous production of ethanol that included cell loading (Taylor et al. in Biotechnol Prog 15:740–751, 2002). The data indicated a minimum production cost at a catalyst loading within 10% of the optimum predicted by the mathematical model. The production rate versus cell loading in most cases displayed a sufficiently broad optimum that the same (non-optimal) rate could be obtained at a significantly smaller loading such as a rate at 80% loading being equal to the rate at 20% loading. These results are particularly important because of the renewed interest in ethanol production (Novozymes and BBI International, Fuel ethanol: a technological evolution, 2004).     

Genome-wide meta-analysis for rheumatoid arthritis

Meta-analysis is being increasingly used as a tool for integrating data from different studies of complex phenotypes, because the power of any one study to identify causal loci is limited. We applied a novel meta-analytical approach (Loesgen et al. in Genet Epidemiol 21(Suppl 1):S142–S147, 2001) in compiling results from four studies of rheumatoid arthritis in Caucasians including two studies from NARAC (Jawaheer et al. in Am J Hum Genet 68:927–936, 2001; Jawaheer et al. in Arthritis Rheum 48:906–916, 2003), one study from the UK (MacKay et al. in Arthritis Rheum 46:632–639, 2001) and one from France (Cornelis et al. in Proc Natl Acad Sci USA 95:10746–10750, 1998). For each study, we obtained NPL scores by performing interval mapping (2 cM intervals) using GeneHunter2 (Kruglyak et al. in Am J Hum Genet 58:1347–1363, 1996; Markianos et al. in Am J Hum Genet 68:963–977, 2001). The marker maps differed among the three consortium groups, therefore, the marker maps were aligned after the interval mapping was completed and the NPL scores that were within 1 cM of each other were combined using the method of Loesgen et al. (Genet Epidemiol 21(Suppl 1):S142–S147, 2001) by calculating the weighted average of the NPL score. This approach avoids some problems in analysis encountered by using GeneHunter2 when some markers in the sample are not genotyped. This procedure provided marginal evidence (P<0.05) of linkage on chromosome 1, 2, 5 and 18, strong evidence (P<0.01) on chromosomes 8 and 16, and overwhelming evidence in the HLA region of chromosome 6.     

Male Mating Interest Varies with Female Fecundity in <Emphasis Type="Italic">Pan troglodytes schweinfurthii</Emphasis> of Kanyawara,Kibale National Park

Female chimpanzees mate promiscuously during a period of extended receptivity marked by prominent sexual swelling. Recent studies of wild chimpanzees indicate that subtle variations in swelling size could act as a reliable cue of female fertilization potential both within and between cycles (Emery and Whitten Behavioral Ecology and Sociobiology, 54, 340–351, 2003; Deschner et al. Hormones and Behavior, 46, 204–215, 2004). Copulation rates increase during the periovulatory period and during conception cycles (Deschner et al. Hormones and Behavior, 46, 204–215, 2004; Emery Thompson American Journal of Primatology, 67, 137–158, 2005a), suggesting that males may be able to assess female fertilization potential. We asked whether facultative timing of copulation in Kanyawara chimpanzees was due to increased male mating interest or to increased female proceptivity during the most fecund days. We assessed multiple measures of male mating effort in cycles aligned relative to the day of detumescence and compared periovulatory days to other days of maximal swelling, and conception cycles to nonconception cycles. The rate and proportion of male initiative in soliciting sexual behavior increased during periods of highest fertilization potential. Males were also more likely to interrupt copulations, associate with estrous females, and compete with other males when females were most likely to conceive. Females initiated copulations more frequently during conception cycles but did not visibly shift mating behavior within cycles. Our results support the hypothesis that male chimpanzees have the ability to assess the profitability of mating attempts, a trait that may act as a counter-adaptation to female strategies to obscure paternity. We discuss potential cues and the implications for female reproductive strategies.     

A general linear framework for the comparison and evaluation of models of sensorimotor synchronization

Sensorimotor synchronization (SMS), the temporal coordination of a rhythmic movement with an external rhythm, has been studied most often in tasks that require tapping along with a metronome. Models of SMS use information about the timing of preceding stimuli and responses to predict when the next response will be made. This article compares the theoretical structure and empirical predictions of four two-parameter models proposed in the literature: Michon (Timing in temporal tracking, Van Gorcum, Assen, 1967), Hary and Moore (Br J Math Stat Psychol 40:109–124, 1987b), Mates (Biol Cybern 70:463–473, 1994a; Biol Cybern 70:475–484, 1994b), and Schulze et al. (Mus Percept 22:461–467, 2005). By embedding these models within a general linear framework, the mathematical equivalence of the Michon, Hary and Moore, and Schulze et al. models is demonstrated. The Mates model, which differs from the other three, is then tested empirically with new data from a tapping experiment in which the metronome alternated between two tempi. The Mates model predictions are found to be invalid for about one-third of the trials, suggesting that at least one of the model’s underlying assumptions is incorrect. The other models cannot be refuted as easily, but they do not predict some features of the data very accurately. Comparison of the models’ predictions in a training/test procedure did not yield any significant differences. The general linear framework introduced here may help in the formulation of new models that make better predictions.