Bayesian Hierarchical Modeling Framework Research Articles

Estimating missing values in China\u2019s official socioeconomic statistics using progressive spatiotemporal Bayesian hierarchical modeling

Due to a large number of missing values, both spatially and temporally, China has not published a complete official socioeconomic statistics dataset at the county level, which is the country’s basic scale of official statistics data collection. We developed a procedure to impute the missing values under the Bayesian hierarchical modeling framework. The procedure incorporates two novelties. First, it takes into account spatial autocorrelations and temporal trends for those easier-to-impute variables with small missing percentages. Second, it further uses the first-step complete variables as covariate information to improve the modeling of more-difficult-to-impute variables with large missing percentages. We applied this progressive spatiotemporal (PST) method to China’s official socioeconomic statistics during 2002–2011 and compared it with four other widely used imputation methods, including k-nearest neighbors (kNN), expectation maximum (EM), singular value decomposition (SVD) and random forest (RF). The results show that the PST method outperforms these methods, thus proving the effects of sophisticatedly incorporating the additional spatial and temporal information and progressively utilizing the covariate information. This study has an outcome that allows China to construct a complete socioeconomic dataset and establishes a methodology that can be generally useful for estimating missing values in large spatiotemporal datasets.

A Bayesian hierarchical model for quantitative and qualitative responses

ABSTRACTIn many science and engineering systems both quantitative and qualitative output observations are collected. If modeled separately the important relationship between the two types of responses is ignored. In this article, we propose a Bayesian hierarchical modeling framework to jointly model a continuous and a binary response. Compared with the existing methods, the Bayesian method overcomes two restrictions. First, it solves the problem in which the model size (specifically, the number of parameters to be estimated) exceeds the number of observations for the continuous response. We use one example to show how such a problem can easily occur if the design of the experiment is not proper; all the frequentist approaches would fail in this case. Second, the Bayesian model can provide statistical inference on the estimated parameters and predictions, whereas it is not clear how to obtain inference using the latest method proposed by Deng and Jin (2015), which jointly models the two responses via constrained likelihood. We also develop a Gibbs sampling scheme to generate accurate estimation and prediction for the Bayesian hierarchical model. Both the simulation and the real case study are shown to illustrate the proposed method.

Risk Assessment and Mapping of Hand, Foot, and Mouth Disease at the County Level in Mainland China Using Spatiotemporal Zero-Inflated Bayesian Hierarchical Models.

Hand, foot, and mouth disease (HFMD) is a worldwide infectious disease, prominent in China. China’s HFMD data are sparse with a large number of observed zeros across locations and over time. However, no previous studies have considered such a zero-inflated problem on HFMD’s spatiotemporal risk analysis and mapping, not to mention for the entire Mainland China at county level. Monthly county-level HFMD cases data combined with related climate and socioeconomic variables were collected. We developed four models, including spatiotemporal Poisson, negative binomial, zero-inflated Poisson (ZIP), and zero-inflated negative binomial (ZINB) models under the Bayesian hierarchical modeling framework to explore disease spatiotemporal patterns. The results showed that the spatiotemporal ZINB model performed best. Both climate and socioeconomic variables were identified as significant risk factors for increasing HFMD incidence. The relative risk (RR) of HFMD at the local scale showed nonlinear temporal trends and was considerably spatially clustered in Mainland China. The first complete county-level spatiotemporal relative risk maps of HFMD were generated by this study. The new findings provide great potential for national county-level HFMD prevention and control, and the improved spatiotemporal zero-inflated model offers new insights for epidemic data with the zero-inflated problem in environmental epidemiology and public health.

Using angling and electric fishing to estimate smallmouth bass abundance in a river

AbstractEstimating abundance is fundamental to effective fishery management but can be challenging in a river where spatial and temporal heterogeneity may preclude the consistent use of a single sampling gear and different gears have differing size selectivity and capture probabilities of fish. In this study, the number of smallmouth bass, Micropterus dolomieu Lacepède, was estimated based on mark–recapture data from angling and boat electric fishing in a 4.2‐km regulated section (mean width = 115 m) of the Broad River, South Carolina, USA. Closed‐population capture–mark–recapture models were fit in the Bayesian hierarchical modelling framework with an estimated number of 2,380 fish (95% credible interval: 1,578–3,693) over 200 mm TL, although simulations indicated that abundance would be slightly overestimated (<20%) when two gears selected for different individuals. Integrating two gear types into a mark–recapture study can provide a method for assessing abundance in spatially or temporally heterogeneous habitats.

Evidence for regional nitrogen stress on chlorophyll a in lakes across large landscape and climate gradients

AbstractNitrogen (N) and phosphorus (P) commonly stimulate phytoplankton production in lakes, but recent observations from lakes from an agricultural region suggest that nitrate may have a subsidy‐stress effect on chlorophyll a (Chl a). It is unclear, however, how generalizable this effect might be. Here, we analyzed a large water quality dataset of 2385 lakes spanning 60 regions across 17 states in the Northeastern and Midwestern U.S. to determine if N subsidy‐stress effects on phytoplankton are common and to identify regional landscape characteristics promoting N stress effects in lakes. We used a Bayesian hierarchical modeling framework to test our hypothesis that Chl a–total N (TN) threshold relationships would be common across the central agricultural region of the U.S. (“the Corn Belt”), where lake N and P concentrations are high. Data aggregated across all regions indicated that high TN concentrations had a negative effect on Chl a in lakes with concurrent high total P. This large‐scale pattern was driven by relationships within only a subset of regions, however. Eight regions were identified as having Chl a–TN threshold relationships, but only two of these regions located within the Corn Belt clearly demonstrated this subsidy‐stress relationship. N stress effects were not consistent across other intense agricultural regions, as we hypothesized. These findings suggest that interactions among regional land use and land cover, climate, and hydrogeology may be important in determining the synergistic conditions leading to N subsidy‐stress effects on lake phytoplankton.

Local overfishing may be avoided by examining parameters of a spatio-temporal model.

Spatial erosion of stock structure through local overfishing can lead to stock collapse because fish often prefer certain locations, and fisheries tend to focus on those locations. Fishery managers are challenged to maintain the integrity of the entire stock and require scientific approaches that provide them with sound advice. Here we propose a Bayesian hierarchical spatio-temporal modelling framework for fish abundance data to estimate key parameters that define spatial stock structure: persistence (similarity of spatial structure over time), connectivity (coherence of temporal pattern over space), and spatial variance (variation across the seascape). The consideration of these spatial parameters in the stock assessment process can help identify the erosion of structure and assist in preventing local overfishing. We use Atlantic cod (Gadus morhua) in eastern Canada as a case study an examine the behaviour of these parameters from the height of the fishery through its collapse. We identify clear signals in parameter behaviour under circumstances of destructive stock erosion as well as for recovery of spatial structure even when combined with a non-recovery in abundance. Further, our model reveals the spatial pattern of areas of high and low density persists over the 41 years of available data and identifies the remnant patches. Models of this sort are crucial to recovery plans if we are to identify and protect remaining sources of recolonization for Atlantic cod. Our method is immediately applicable to other exploited species.

Detecting long-term occupancy changes in Californian odonates from natural history and citizen science records

In a world of rapid environmental change, effective biodiversity conservation and management relies on our ability to detect changes in species occurrence. While long-term, standardized monitoring is ideal for detecting change, such monitoring is costly and rare. An alternative approach is to use historical records from natural history collections as a baseline to compare with recent observations. Here, we combine natural history collection data with citizen science observations within a hierarchical Bayesian occupancy modeling framework to identify changes in the occupancy of Californian dragonflies and damselflies (Odonata) over the past century. We model changes in the probability of occupancy of 34 odonate species across years and as a function of climate, after correcting for likely variation in detection probability using proxies for recorder effort and seasonal variation. We then examine whether biological traits can help explain variation in temporal trends. Models built using only opportunistic records identify significant changes in occupancy across years for 14 species, with eight of those showing significant declines and six showing significant increases in occupancy in the period 1900–2013. These changes are consistent with estimates obtained using more standardized resurvey data, regardless of whether resurvey data are used individually or in conjunction with the opportunistic dataset. We find that species increasing in occupancy over time are also those whose occupancy tends to increase with higher minimum temperatures, which suggests that these species may be benefiting from increasing temperatures across California. Furthermore, these species are also mostly habitat generalists, whilst a number of habitat specialists display some of the largest declines in occupancy across years. Our approach enables more robust estimates of temporal trends from opportunistic specimen and observation data, thus facilitating the use of these data in biodiversity conservation and management.

Data Integration Model for Air Quality: A Hierarchical Approach to the Global Estimation of Exposures to Ambient Air Pollution

SummaryAir pollution is a major risk factor for global health, with 3 million deaths annually being attributed to fine particulate matter ambient pollution (PM2.5). The primary source of information for estimating population exposures to air pollution has been measurements from ground monitoring networks but, although coverage is increasing, regions remain in which monitoring is limited. The data integration model for air quality supplements ground monitoring data with information from other sources, such as satellite retrievals of aerosol optical depth and chemical transport models. Set within a Bayesian hierarchical modelling framework, the model allows spatially varying relationships between ground measurements and other factors that estimate air quality. The model is used to estimate exposures, together with associated measures of uncertainty, on a high resolution grid covering the entire world from which it is estimated that 92% of the world's population reside in areas exceeding the World Health Organization's air quality guidelines.

Socio-environmental and measurement factors drive variation in influenza-like illness

ObjectiveTo assess the use of medical claims records for surveillance andepidemiological inference through a case study that examines howecological and social determinants and measurement error contributeto spatial heterogeneity in reports of influenza-like illness across theUnited States.IntroductionTraditional infectious disease epidemiology is built on thefoundation of high quality and high accuracy data on disease andbehavior. Digital infectious disease epidemiology, on the other hand,uses existing digital traces, re-purposing them to identify patterns inhealth-related processes. Medical claims are an emerging digital datasource in surveillance; they capture patient-level data across an entirepopulation of healthcare seekers, and have the benefits of medicalaccuracy through physician diagnoses, and fine spatial and temporalresolution in near real-time.Our work harnesses the large volume and high specificity ofdiagnosis codes in medical claims to improve our understanding ofthe mechanisms driving spatial variation in reported influenza activityeach year. The mechanisms hypothesized to drive these patterns areas varied as: environmental factors affecting transmission or virussurvival, travel flows between different populations, population agestructure, and socioeconomic factors linked to healthcare access andquality of life. Beyond process mechanisms, the nature of surveillancedata collection may affect our interpretation of spatial epidemiologicalpatterns [1], particularly since influenza is a non-reportable diseasewith non-specific symptoms ranging from asymptomatic to severe.Considering the ways in which medical claims are generated, biasesmay arise from healthcare-seeking behavior, insurance coverage, andmedical claims database coverage in study populations.MethodsUsing aggregated U.S. medical claims for influenza-like illness(ILI) from the 2001-2002 through 2008-2009 flu seasons [2],we developed a Bayesian hierarchical modeling framework toestimate the importance of both ecological and social determinantsand measurement-related factors on observed county-level variationof influenza disease burden across the United States. Integrated NestedLaplace Approximation (INLA) techniques for Bayesian inferencewere used to render our questions computationally tractable due tothe high spatial resolution of our data (Figure 1) and the multiplicityof models in our analysis [3]. Linking data from a variety of publiclyavailable sources, we determined the strength, directionality, andconsistency of these factors over multiple flu seasons.ResultsWe found that measurement-related factors – healthcare-seekingbehavior, insurance coverage, and medical claims database coverage– were strong predictors of greater ILI intensity across seasons.Secondarily, poverty and specific humidity were negatively associatedwith ILI intensity for several seasons. Finally, by incorporatingmechanistic and measurement factors into our model, our modelpredictions present an improved map of influenza-like illness in theUnited States for the flu seasons in our study period.ConclusionsWe present a flexible modeling approach that applies to differentmedical claims diagnosis codes and disease surveillance data anddemonstrates the utility of Bayesian hierarchical models for large-scale ecological analyses. Our results increase our knowledge of thespatial distribution of influenza and the underlying processes thatdrive these patterns, promote finer spatial targeting for differenttypes of interventions, and enable the interpolation of burden in areasdifficult to surveil through traditional public health. Moreover, theyhighlight the relative contributions of surveillance data collectionand ecological processes to spatial variation in disease, and highlightthe importance of considering measurement biases when usingsurveillance data for epidemiological inference.

A Bayesian hierarchical model for estimating spatial and temporal variation in vegetation phenology from Landsat time series

Phenology is a key indicator of vegetation response to global change. Satellite observations from coarse resolution sensors (i.e., MODIS and AVHRR) have been widely used to study impacts of climate anomalies on vegetation phenology. The advantage of coarse resolution sensors are daily observations across the entire globe, but the coarse spatial resolution, as well as the relatively short time span covered by MODIS, are significant drawbacks in analyzing landscape-scale trends in phenology. Time series data from the medium resolution sensor Landsat may overcome these issues. However, because of Landsat's lower observation frequency, phenological methods developed with coarse resolution data are not directly transferable. Here, we demonstrate a new Bayesian hierarchical modeling framework for estimating inter-annual variation in vegetation phenology from Landsat time series while controlling for spatial variation. The method pools all available observations to estimate the spatial variation in phenological parameters, while specifically modeling inter-annual variation as random effect terms. The advantage of a Bayesian approach is the ability to incorporate prior knowledge from other phenology and climate observations to reduce variability in the estimates, as well as a more robust estimation of uncertainty. We demonstrate and evaluate the modeling framework with a case study of changing spring phenology in broad-leaved trees in the Bavarian Forest National Park in southern Germany. Results show that the model estimated the spatial and temporal variation in phenological parameters precisely. Temporal variation in start of season showed overall strong agreement with ground-based measures of bud-break variability (r=0.82 [0.80–0.84]). Our proposed modeling framework will help to better monitor and understand changes in vegetation phenology at scales yet unexplored by the phenological community.

A meta-analysis of fecundity in rockfishes (genus Sebastes)

Estimates of fecundity at size are a key component of age-structured population dynamics models such as those used for stock assessment purposes. Fecundity-length relationships for many exploited fish stocks are often poorly understood due to small sample sizes and high levels of intrinsic variability among individuals. Several species in the genus Sebastes match this description, and are important components of commercial and recreational fisheries off the U.S. West Coast. We used a hierarchical Bayesian modeling framework for simultaneous estimation of parameters in the fecundity-length relationships for 29 species of rockfish while accounting for variability within and among subgenera. Hierarchical models that account for subgeneric phylogenetic structure and variability among species outperformed models that treat all species in the genus as exchangeable or ignored variation among species within a subgenus, based on a predictive information criterion. Our results confirm that weight-specific fecundity increases with size in nearly all Sebastes species, and that the assumption of proportionality between mature female biomass and total egg production is inappropriate for most Sebastes stock assessments. Benefits of the hierarchical framework include prediction of fecundity-length parameters at the genus, subgenus, and species level, as well as predictive distributions for unobserved Sebastes species.

Modeling Day-to-Day Variability of Glucose-Insulin Regulation Over 12-Week Home Use of Closed-Loop Insulin Delivery.

Parameters of physiological models of glucose-insulin regulation in type 1 diabetes have previously been estimated using data collected over short periods of time and lack the quantification of day-to-day variability. We developed a new hierarchical model to relate subcutaneous insulin delivery and carbohydrate intake to continuous glucose monitoring over 12 weeks while describing day-to-day variability. Sensor glucose data sampled every 10-min, insulin aspart delivery and meal intake were analyzed from eight adults with type 1 diabetes (male/female 5/3, age 39.9 ± 9.5 years, BMI 25.4 ± 4.4kg/m2, HbA1c 8.4 ± 0.6% ) who underwent a 12-week home study of closed-loop insulin delivery. A compartment model comprised of five linear differential equations; model parameters were estimated using the Markov chain Monte Carlo approach within a hierarchical Bayesian model framework. Physiologically, plausible a posteriori distributions of model parameters including insulin sensitivity, time-to-peak insulin action, time-to-peak gut absorption, and carbohydrate bioavailability, and good model fit were observed. Day-to-day variability of model parameters was estimated in the range of 38-79% for insulin sensitivity and 27-48% for time-to-peak of insulin action. In conclusion, a linear Bayesian hierarchical approach is feasible to describe a 12-week glucose-insulin relationship using conventional clinical data. The model may facilitate in silico testing to aid the development of closed-loop insulin delivery systems.

Spatio-temporal assimilation of modelled catchment loads with monitoring data in the Great Barrier Reef

Soil erosion and sediment transport into waterways and the ocean can adversely affect water clarity, leading to the deterioration of marine ecosystems such as the iconic Great Barrier Reef (GBR) in Australia. Quantifying a sediment load and its associated uncertainty is an important task in delineating how changes in management practices can contribute to improvements in water quality, and therefore continued sustainability of the GBR. However, monitoring data are spatially (and often temporally) sparse, making load estimation complicated, particularly when there are lengthy periods between sampling or during peak flow periods of major events when samples cannot be safely taken. We develop a spatio-temporal statistical model that is mechanistically motivated by a process-based deterministic model called Dynamic SedNet. The model is developed within a Bayesian hierarchical modelling framework that uses dimension reduction to accommodate seasonal and spatial patterns to assimilate monitored sediment concentration and flow data with output from Dynamic SedNet. The approach is applied in the Upper Burdekin catchment in Queensland, Australia, where we obtain daily estimates of sediment concentrations, stream discharge volumes and sediment loads at 411 spatial locations across 20 years. Our approach provides a method for assimilating both monitoring data and modelled output, providing a statistically rigorous method for quantifying uncertainty through space and time that was previously unavailable through process-based models.

Integrating independent spatio-temporal replications to assess population trends in disease spread.

Glaucoma is the second leading cause of blindness in the USA. A visual field test (perimetry) is used to sample and quantitate visual field function in preselected regions in the eye. These regions can be considered a spatial field with replications across independently measured individuals. At return visits, a new set of visual field measurements is obtained producing a subject specific spatio-temporal dataset. We develop a Bayesian hierarchical modeling framework to analyze these spatio-temporal datasets both for individual level spread and as aggregate population level trends. Our model extends previous research utilizing a dimension reduction matrix and individual specific latent variables. Human characteristics are incorporated into the model to help explain glaucoma progression. One beneficial product of our model is smoothed estimates for individuals. We also specify how progression rates are computed for monitoring purposes so that clinicians can track changes and predict forward in time. Copyright © 2016 John Wiley & Sons, Ltd.

Combining growth curves when a longitudinal study switches measurement tools.

When longitudinal studies are performed to investigate the growth of traits in children, the measurement tool being used to quantify the trait may need to change as the subjects' age throughout the study. Changing the measurement tool at some point in the longitudinal study makes the analysis of that growth challenging which, in turn, makes it difficult to determine what other factors influence the growth rate. We developed a Bayesian hierarchical modeling framework that relates the growth curves per individual for each of the different measurement tools and allows for covariates to influence the shapes of the curves by borrowing strength across curves. The method is motivated by and demonstrated by speech perception outcome measurements of children who were implanted with cochlear implants. Researchers are interested in assessing the impact of age at implantation and comparing the growth rates of children who are implanted under the age of two versus those implanted between the ages of two and four.

Is the cluster environment quenching the Seyfert activity in elliptical and spiral galaxies?

We developed a hierarchical Bayesian model (HBM) to investigate how the presence of Seyfert activity relates to their environment, herein represented by the galaxy cluster mass, $M_{200}$, and the normalized cluster-centric distance, $r/r_{200}$. We achieved this by constructing an unbiased sample of galaxies from the Sloan Digital Sky Survey, with morphological classifications provided by the Galaxy Zoo Project. A propensity score matching approach is introduced to control for the effects of confounding variables: stellar mass, galaxy colour, and star formation rate. The connection between Seyfert-activity and environmental properties in the de-biased sample is modelled within an HBM framework using the so-called logistic regression technique, suitable for the analysis of binary data (e.g., whether or not a galaxy hosts an AGN). Unlike standard ordinary least square fitting methods, our methodology naturally allows modelling the probability of Seyfert-AGN activity in galaxies on their natural scale, i.e. as a binary variable. Furthermore, we demonstrate how an HBM can incorporate information of each particular galaxy morphological type in a unified framework. In elliptical galaxies, our analysis indicates a strong correlation of Seyfert-AGN activity with $r/r_{200}$, and a weaker correlation with the mass of the host. In spiral galaxies these trends do not appear, suggesting that the link between Seyfert activity and the properties of spiral galaxies are independent of the environment.

Bayesian spatially dependent variable selection for small area health modeling.

Statistical methods for spatial health data to identify the significant covariates associated with the health outcomes are of critical importance. Most studies have developed variable selection approaches in which the covariates included appear within the spatial domain and their effects are fixed across space. However, the impact of covariates on health outcomes may change across space and ignoring this behavior in spatial epidemiology may cause the wrong interpretation of the relations. Thus, the development of a statistical framework for spatial variable selection is important to allow for the estimation of the space-varying patterns of covariate effects as well as the early detection of disease over space. In this paper, we develop flexible spatial variable selection approaches to find the spatially-varying subsets of covariates with significant effects. A Bayesian hierarchical latent model framework is applied to account for spatially-varying covariate effects. We present a simulation example to examine the performance of the proposed models with the competing models. We apply our models to a county-level low birth weight incidence dataset in Georgia.

Fetal growth restriction and the immunological consequences of mother re-exposure to orally tolerated proteins

With more than a million spectators expected to travel among 12 different cities in Brazil during the football World Cup, June 12–July 13, 2014, the risk of the mosquito-transmitted disease dengue fever is a concern. We addressed the potential for a dengue epidemic during the tournament, using a probabilistic forecast of dengue risk for the 553 microregions of Brazil, with risk level warnings for the 12 cities where matches will be played.We obtained real-time seasonal climate forecasts from several international sources (European Centre for Medium-Range Weather Forecasts [ECMWF], Met Office, Meteo-France and Centro de Previsão de Tempo e Estudos Climáticos [CPTEC]) and the observed dengue epidemiological situation in Brazil at the forecast issue date as provided by the Ministry of Health. Using this information we devised a spatiotemporal hierarchical Bayesian modelling framework that enabled dengue warnings to be made 3 months ahead. By assessing the past performance of the forecasting system using observed dengue incidence rates for June, 2000–2013, we identified optimum trigger alert thresholds for scenarios of medium-risk and high-risk of dengue.Our forecasts for June, 2014, showed that dengue risk was likely to be low in the host cities Brasília, Cuiabá, Curitiba, Porto Alegre, and São Paulo. The risk was medium in Rio de Janeiro, Belo Horizonte, Salvador, and Manaus. High-risk alerts were triggered for the northeastern cities of Recife (phigh=19%), Fortaleza (phigh=46%), and Natal (phigh=48%). For these high-risk areas, particularly Natal, the forecasting system did well for previous years (in June, 2000–13).This timely dengue early warning permits the Ministry of Health and local authorities to implement appropriate, city-specific mitigation and control actions ahead of the World Cup.European Commission's Seventh Framework Research Programme projects DENFREE, EUPORIAS, and SPECS; Conselho Nacional de Desenvolvimento Científico e Tecnológico and Fundação de Amparo à Pesquisa do Estado do Rio de Janeiro.

Interception of wet deposited atmospheric pollutants by herbaceous vegetation: Data review and modelling

Better understanding and predicting interception of wet deposited pollutants by vegetation remains a key issue in risk assessment studies of atmospheric pollution. We develop different alternative models, following either empirical or semi-mechanistic descriptions, on the basis of an exhaustive dataset consisting of 440 observations obtained in controlled experiments, from 1970 to 2014, for a wide variety of herbaceous plants, radioactive substances and rainfall conditions. The predictive performances of the models and the uncertainty/variability of the parameters are evaluated under Hierarchical Bayesian modelling framework. It is demonstrated that the variability of the interception fraction is satisfactorily explained and quite accurately modelled by a process-based alternative in which absorption of ionic substances onto the foliage surfaces is determined by their electrical valence. Under this assumption, the 95% credible interval of the predicted interception fraction encompasses 81% of the observations, including situations where either plant biomass or rainfall intensity are unknown. This novel approach is a serious candidate to challenge existing empirical relationships in radiological or chemical risk assessment tools.

Bayesian regularized regression based on composite quantile method

Recently, variable selection based on penalized regression methods has received a great deal of attention, mostly through frequentist’s models. This paper investigates regularization regression from Bayesian perspective. Our new method extends the Bayesian Lasso regression (Park and Casella, 2008) through replacing the least square loss and Lasso penalty by composite quantile loss function and adaptive Lasso penalty, which allows different penalization parameters for different regression coefficients. Based on the Bayesian hierarchical model framework, an efficient Gibbs sampler is derived to simulate the parameters from posterior distributions. Furthermore, we study the Bayesian composite quantile regression with adaptive group Lasso penalty. The distinguishing characteristic of the newly proposed method is completely data adaptive without requiring prior knowledge of the error distribution. Extensive simulations and two real data examples are used to examine the good performance of the proposed method. All results confirm that our novel method has both robustness and high efficiency and often outperforms other approaches.