Hamiltonian Monte Carlo method for estimating variance components

A Hamiltonian Monte Carlo algorithm is a Markov chain Monte Carlo method, and the method has a potential to improve estimating parameters effectively. Hamiltonian Monte Carlo is based on Hamiltonian dynamics, and it follows Hamilton's equations, which are expressed as two differential equations. In the sampling process of Hamiltonian Monte Carlo, a numerical integration method called leapfrog integration is used to approximately solve Hamilton's equations, and the integration is required to set the number of discrete time steps and the integration stepsize. These two parameters require some amount of tuning and calibration for effective sampling. In this study, we applied the Hamiltonian Monte Carlo method to animal breeding data and identified the optimal tunings of leapfrog integration for normal and inverse chi-square distributions. Then, using real pig data, we revealed the properties of the Hamiltonian Monte Carlo method with the optimal tuning by applying models including variance explained by pedigree information or genomic information. Compared with the Gibbs sampling method, the Hamiltonian Monte Carlo method had superior performance in both models. We have provided the source codes of this method written in the Fortran language at https://github.com/A-ARAKAWA/HMC .     

Parametric bootstrap for testing model fitting in the proportional hazards framework: an application to the survival analysis of Bruna dels Pirineus beef calves

Given that correct assumptions on the baseline survival function are determinant for the validity of further inferences, specific tools to test the fit of a model to real data become essential in proportional hazards models. In this sense, we have proposed a parametric bootstrap to test the fit of survival models. Monte Carlo simulations are used to generate new data sets from the estimates obtained through the assumed models, and then bootstrap intervals can be established for the survival function along the time space studied. Significant fitting deficiencies are revealed when the real survival function is not included within the bootstrap interval. We tested this procedure in a survival data set of Bruna dels Pirineus beef calves, assuming 4 parametric models (exponential, Weibull, exponential time-dependent, Weibull time-dependent) and the Cox's semiparametric model. Fitting deficiencies were not observed for the Cox's model and the exponential time-dependent model, whereas the Weibull time-dependent model suffered from moderate overestimation at different ages. Thus, the exponential time-dependent model appears to be preferable because of its correct fit for survival data of beef calves and its smaller computational and time requirements. Exponential and Weibull models were completely rejected due to the continuous over- and underestimation of the survival probability reported. Results here highlighted the flexibility of parametric models with time-dependent effects, achieving a fit comparable to nonparametric models.     

Implementation issues for Markov Chain Monte Carlo methods in random regression test-day models

Markov Chain Monte Carlo methods made possible estimation of parameters for complex random regression test‐day models. Models evolved from single‐trait with one set of random regressions to multiple‐trait applications with several random effects described by regressions. Gibbs sampling has been used for models with linear (with respect to coefficients) regressions and normality assumptions for random effects. Difficulties associated with implementations of Markov Chain Monte Carlo schemes include lack of good practical methods to assess convergence, slow mixing caused by high posterior correlations of parameters and long running time to generate enough posterior samples. Those problems are illustrated through comparison of Gibbs sampling schemes for single‐trait random regression test‐day models with different model parameterizations, different functions used for regressions and posterior chains of different sizes. Orthogonal polynomials showed better convergence and mixing properties in comparison with ‘lactation curve’ functions of the same number of parameters. Increasing the order of polynomials resulted in smaller number of independent samples for covariance components. Gibbs sampling under hierarchical model parameterization had a lower level of autocorrelation and required less time for computation. Posterior means and standard deviations of genetic parameters were very similar for chains of different size (from 20 000 to 1 000 000) after convergence. Single‐trait random regression models with large data sets can be analysed by Markov Chain Monte Carlo methods in relatively short time. Multiple‐trait (lactation) models are computationally more demanding and better algorithms are required.     

The use of multilevel models to evaluate sources of variation in reproductive performance in dairy cattle in Reunion Island.

Sources of variation in measures of reproductive performance in dairy cattle were evaluated using data collected from 3207 lactations in 1570 cows in 50 herds from five geographic regions of Reunion Island (located off the east coast of Madagascar). Three continuously distributed reproductive parameters (intervals from calving-to-conception, calving-to-first-service and first-service-to-conception) were considered, along with one Binomial outcome (first-service-conception risk). Multilevel models which take into account the hierarchical nature of the data were used to fit all models. For the overall measure of calving-to-conception interval, 86% of the variation resided at the lactation level with only 7, 6 and 2% at the cow, herd and regional levels, respectively. The proportion of variance at the herd and cow levels were slightly higher for the calving-to-first-service interval (12 and 9%, respectively) - but for the other two parameters (first-service-conception risk and first-service-to-conception interval), >90% of the variation resided at the lactation level. For the three continuous dependent variables, comparison of results between models based on log-transformed data and Box-Cox-transformed data suggested that minor departures from the assumption of normality did not have a substantial effect on the variance estimates. For the Binomial dependent variable, five different estimation procedures (penalised quasi-likelihood, Markov-Chain Monte Carlo, parametric and non-parametric bootstrap estimates and maximum-likelihood) yielded substantially different results for the estimate of the cow-level variance.     

The DairyCHAMP program: a computerised recording system for dairy herds

The DairyCHAMP program is an animal health and management software program that helps daily animal management, herd performance monitoring and problem analysis. Data entry to the program uses a data dictionary and includes an error-checking system that ensures the consistency and appropriateness of data entered. DairyCHAMP performs health management functions, provides a convenient user interface, ensures uniform data across farms by using a standard data dictionary, can be fully integrated with decision-making software programs like DairyORACLE, and is flexible enough to be useful for many types of dairy facilities. Data are entered via a menu-based system. Animal events are organised around reproduction and lactation cycles and health records. Farm records include inventories for drug, feed and semen. Farm parameters can be established which customize the program for an individual farm. The database system is an integration of three schemas: the individual user's view, the community view and the storage system. The individual user's view must be easy to use, while the storage system must be compact enough to fit within the disc storage space on a microcomputer. This conflict requires a translation from one schema to another. The DairyCHAMP program accomplishes this through a coding system which assigns a code number to each event. The program can add synonyms to this event dictionary by assigning the same code number to the synonym the user chooses. The DairyCHAMP program provides access to the large amounts of data required to aid in daily animal management, allow performance monitoring and analyse problems. Its highly integrated system is efficient and easy to use and maintain.     

Seasonal variation of bulk milk somatic cell counts in UK dairy herds: investigations of the summer rise

Individual cow somatic cell count (SCC) patterns were explored over a one year period in 33 dairy herds to investigate the reason for a summer rise in bulk milk somatic cell counts (BMSCC). Cow test day somatic cell counts were categorised according to the magnitude of change since the previous test day reading, to examine which categories were responsible for the summer increase. Multilevel models using Markov chain Monte Carlo methods were specified to estimate the number of somatic cells/ml produced by different cell count categories. Stage of lactation and parity were accounted for in the models. There was an increase in the proportion of cows that remained above 200,000 cells/ml for two consecutive recordings in summer and this group of cows were responsible for 70.8% of the increase in somatic cells/ml produced from May to September compared with October to March. There was no evidence that a greater new infection rate (somatic cell counts moving from below 100,000 cells/ml to over 200,000 cells/ml) contributed to the increased summer bulk milk somatic cell counts. There was no indication that a general small increase in all somatic cell counts played an important role in the increased summer somatic cell counts. Markov chain Monte Carlo methods provided a valuable and flexible platform for parameter estimation in reasonably complex multilevel models.     

Monte Carlo模拟方法在动物育种中的应用

Monte Carlo方法是一类通过随机变量统计试验,随机模拟以求得问题近似解的方法。已有的研究结果表明,利用计算机模拟方法研究动物遗传育种理论与实践问题不仅是可行的,而且是可靠的。应用计算机模拟动物育种问题时需要建立一个准确可靠的遗传传递模型,利用能产生周期长且统计性质优的伪随机数产生方法,以使模拟结果更加准确可靠。它在动物育种方案研究、人工选择、QTL定位及保种方案等方面的模拟结果,为育种理论研究和实际育种工作提供了强有力的支撑。     

Estimating danofloxacin withdrawal time in broiler chickens based on physiologically based pharmacokinetics modeling

In this study, a physiologically based pharmacokinetics (PBPK) model was firstly developed for danofloxacin in healthy broiler chickens after a single oral administration at 5 mg/kg bw. Then, the model extrapolation from healthy chickens to those infected with Pasteurella multocidaones was performed. The healthy model was validated through a comparison of predicted and previously published concentrations, which indicated that the healthy PBPK model had good predictive ability in plasma, lung, muscle, liver, and kidney, especially at the later sampling time points. Multiple dosing of administration was incorporated into the healthy and infected models. In addition, a Monte Carlo simulation (MCS) included 1000 iterations was further incorporated into both models to predict the withdrawal times of danofloxacin in healthy and infected chickens, which were estimated to be 3 and 2 days, respectively.     

Evaluation of a contract breeding management program in selected Ohio dairy herds with event-time analysis II. Parametric frailty models

The effect of a contract breeding program offered by a breeding co-operative was assessed using parametric frailty models with event-time analysis technique in a field study of Ohio dairies. The program featured tail chalking and daily evaluation of cows for insemination by co-operative technicians; dairy employees no longer handled estrus detection activities. Test day records were obtained between early 2002 and mid-2004 for 16,453 lactations representing 11,398 cows in 31 herds identified as well-managed client herds by the breeding co-operative. Various parametric distributions for event times available in a commercial software (Stata 9.1, College Station, TX) were tested to assess which distribution fit the calving-to-conception data best. After identifying the distribution with the best fit, a full model with potential confounders and other significant predictors of time to pregnancy was developed and then frailty terms were included in the model. Generalized gamma and log-normal distributions fit the data best, but since gamma distribution does not allow the use of frailty effects, log-normal distribution was used in further modeling. Separate accelerated failure time models with frailty terms to account for latent effects at the herd, cow, or lactation level were developed, testing both gamma and inverse Gaussian frailty distributions. In these models, potential confounders and statistically significant predictors were also controlled for, and the association between the contract breeding program and the mean time to pregnancy was characterized using time ratios. The log-normal model identified that interval to pregnancy was associated with breed, herd size, use of ovulation synchronization protocols, parity, calving season and somatic cell score (above or below 4.5) and maximum milk yield prior to pregnancy or censoring. While controlling for these factors, there was a reduction in average time to pregnancy among cows managed under the contract breeding program. All frailty terms were highly significant, regardless of whether it was an individual frailty at the lactation level or a shared frailty at the cow or herd level, suggesting that there was considerable heterogeneity within these levels. Inclusion of a frailty term at the herd level changed the estimate for the contract breeding program considerably, while a frailty term on other levels did not, indicating that herd characteristics (e.g., overall management) have a substantial impact on reproductive performance and should be accounted for in the analysis. Interpretation using time ratios with or without a shared herd frailty found that the contract breeding program was associated with a reduction of 6.5% and 14.1% in mean time to pregnancy, respectively.     

A simulated surveillance program for bovine paratuberculosis in dairy herds in Norway

Monte Carlo simulation models were used to evaluate the feasibility and potential results of a proposed national survey of the prevalence of bovine paratuberculosis (PTB) in dairy herds in Norway. The expected herd prevalence was assumed to be 0.2% in the simulations. The low sensitivity of the ELISA test, the assumed low herd prevalence, the typical low within-herd prevalence of PTB and the small herd sizes all present problems in detection of the disease. Simulations with 500, 1000, 2500 and 6000 herds tested were done. Our results suggest that a national survey would not be feasible at present, due to the low probability of detecting infected herds and because of the high number of false-positive reactions that would be expected to occur.     

Genetic analysis of 'breeding field test status' in Icelandic horses

Genetic evaluation of Icelandic horses is currently based on results from breeding field tests where riding ability and conformation of the horses are evaluated over the course of 1-2 days. Only a small part of registered horses attend these field tests, and it can be assumed that these are not a random sample of the population. In this study, the trait test status was introduced, describing whether a horse was assessed in a breeding field test. This trait was analysed to find out whether it has a genetic variation and how it correlates genetically to other traits in the breeding goal. Breeding field test data included 39,443 mares born in Iceland in 1990-2001, of which 7431 were assessed in the period 1994-2007. The trait was defined in relation to age, gender and stud of horses. Variance and covariance components were estimated using the Markov Chain Monte Carlo method by applying the Gibbs sampler procedure in the DMU program. Three multivariate analyses were performed where the test status trait was analysed with breeding field test traits. Animal models and sire models were applied. Based on estimated heritabilities (0.51-0.67) and genetic correlations (0.00-0.87), the test status trait showed significant genetic variation and was strongly correlated to some traits. The test status trait reflects preselection in the breeding field test traits and should be included in the genetic evaluation to enhance the procedure, reduce selection bias and increase accuracy of the estimation.     

A retrospective analysis of the infectious bovine rhinotracheitis (bovine herpes virus-1) surveillance program in Norway using Monte Carlo simulation models.

Serological surveillance for antibodies against bovine herpes virus type I (BHV-1) which causes infectious bovine rhinotracheitis and infectious pustular vulvovaginitis has been carried out since 1992 in Norway. Since 1993 (when a single infected herd was detected) all bulk-milk and pooled-serum samples have been negative for BHV-1 antibodies. This paper describes the use of Monte Carlo simulation models for the analysis and interpretation of the results of the surveillance and provides support for the contention that the Norwegian cattle population is not infected by BHV-1.     

Community people-pet programs that work

Excellent models exist for people-pet programs in institutions and in the community. Veterinarians should assess the needs of their local communities and adapt a model program to fit these needs.     

Development of a monitoring program for antibiotic resistance in bacteria from Swiss food-producing animals

Resistant bacteria from food-producing animals may compromise the success of antibiotic treatment in animals and in humans. Therefore, the level of antibiotic resistance in bacteria from farm animals and its development over time needs to be monitored. In Switzerland, a monitoring program for antibiotic resistance is currently being developed. Pilot-monitoring programs were conducted in selected animal species in order to obtain current data on antibiotic resistance. The data on the prevalence of bacteria and antibiotic resistance in poultry were used to optimize the sampling plan. The influence of sampling more farms compared to sampling more animals per farm on the prevalence estimate for antibiotic resistance was analyzed by a Monte Carlo simulation model. Accounting for the costs for sample collection, transportation and laboratory analysis of the samples, the number of samples to be taken at the respective step in the production line was optimized. Optimization was defined as maximizing the precision of the prevalence estimate while minimizing the costs. The model will be expanded to other bacterial and animal species in the future.     

Bayesian inference about parameters of a longitudinal trajectory when selection operates on a correlated trait

A hierarchical model for inferring the parameters of the joint distribution of a trait measured longitudinally and another assessed cross-sectionally, when selection has been applied to the cross-sectional trait, is presented. Distributions and methods for a Bayesian implementation via Markov Chain Monte Carlo procedures are discussed for the case where information about the selection criterion is available for all the individuals, but longitudinal records are available only in the later generations. Alternative specifications of the residual covariance structure are suggested. The procedure is illustrated with an analysis of correlated responses in growth curve parameters in a population of rabbits selected for increased growth rate. Results agree with those obtained in a previous study using both selected and control populations. The high correlation between samples indicates slow mixing, resulting in small effective sample sizes and high Monte Carlo standard errors.     

Forecasting clinical disease in pigs: comparing a naive and a Bayesian approach

In veterinary practice the clinician often evaluates and predicts herd health status over time according to clinical criteria. In this paper, we modeled three different clinical signs among pigs based on longitudinal clinical observations in 15 pig herds. We compared and discussed the outputs from two different approaches for making clinical forecasts in a herd: a naive approach using a simple time series model with previous disease observations as predictors and a Bayesian state space models approach, in which the time lag variable entered into the random component of the model. We used the Markov chain Monte Carlo technique to calculate posterior distributions of the forecasts. For the herd specific forecasts the results showed that there were only minor differences between the forecasts from the simple time series model and the median forecasts from the Bayesian model. However, the credibility intervals from the Bayesian model were wider than the forecasts from the simple model and, therefore the Bayesian model encompassed the variability in the forecasts better. Compared to the statistical model, the simple time series would be easier to implement in a practical setting. However, the latter lacks the inherent “generality” from the statistical model that allows the user to make statements about the distribution of the herds and to predict disease status based on the “average” correlation among the herds. The applicability of the Bayesian approach within a clinical decision-making framework was discussed, with special emphasis on the use of prior information and clinical forecasting.     

蒙特卡罗模拟法在抗微生物药物药动学和药效学研究中的应用

常见病原微生物对抗微生物药物的耐药性正逐渐增加,为了达到最佳治疗效果,临床用药必须根据药动学与药效学数据调整给药方案。药动学能够提供药物浓度在组织、体液和感染部位的经时过程,而药效学则反映药物对致病菌的杀灭或抑制能力。蒙特卡罗模拟法则是利用统计学抽样来获得数学方程的近似解的一种方法,目前采用蒙特卡罗模拟法进行实时模拟正成为国际上研究抗微生物药物的药动学和药效学的热点。论文就蒙特卡罗模拟法的原理、拟合过程及其在估算细菌对药物的敏感性折点、比较药动-药效参数以选择最优药物等方面做一综述。     

Simulating the effect of eliminating routine bacterial screening on the negative predictive value of the Ontario fish hatchery disease monitoring program.

We determined the impact of eliminating routine screening for Aeromonas salmonicida and Yersinia ruckeri on the efficacy of the Ontario Ministry of Natural Resources (OMNR) fish disease monitoring program, using Monte Carlo simulation. Because the main purpose of the program is to prevent transferring infected fish among OMNR hatcheries, or to wild fish populations through stocking waterways, the hatchery-level negative predictive value (HNPV) was used as an indicator of monitoring efficacy. The present program (which includes both routine screening of asymptomatic hatchery fish, and diagnostic testing of hatchery mortalities and clinically diseased fish) was confirmed to have a high median HNPV (0.999) for both study pathogens. Simulations suggested that the median probabilities that a hatchery would be pathogen-free if only diagnostic testing were continued (i.e. if no asymptomatic lots were screened), and all diseased lots tested negative for A. salmonicida and Y. ruckeri would be 0.994 for both pathogens (with <5% probability that HNPV would be less than 0.953 and 0.957, respectively) - indicating acceptable monitoring efficacy. However, limitations of the theoretical monitoring model must be considered.     

Pharmacokinetics,milk penetration and PK/PD analysis by Monte Carlo simulation of marbofloxacin,after intravenous and intramuscular administration to lactating goats

The main objectives of this study were (i) to evaluate the serum pharmacokinetic behaviour and milk penetration of marbofloxacin (MFX; 5 mg/kg), after intravenous (IV) and intramuscular (IM) administration in lactating goats and simulate a multidose regimen on steady‐state conditions, (ii) to determine the minimum inhibitory concentration (MIC) and mutant prevention concentration (MPC) of coagulase negative staphylococci (CNS) isolated from caprine mastitis in Córdoba, Argentina and (iii) to make a PK/PD analysis by Monte Carlo simulation from steady‐state pharmacokinetic parameters of MFX by IV and IM routes to evaluate the efficacy and risk of the emergence of resistance. The study was carried out with six healthy, female, adult Anglo Nubian lactating goats. Marbofloxacin was administered at 5 mg/kg bw by IV and IM route. Serum and milk concentrations of MFX were determined with HPLC/uv. From 106 regional strains of CNS isolated from caprine mastitis in herds from Córdoba, Argentina, MICs and MPCs were determined. MIC₉₀ and MPC₉₀ were 0.4 and 6.4 μg/ml, respectively. MIC and MPC‐based PK/PD analysis by Monte Carlo simulation indicates that IV and IM administration of MFX in lactating goats may not be adequate to recommend it as an empirical therapy against CNS, because the most exigent endpoints were not reached. Moreover, this dose regimen could increase the probability of selecting mutants and resulting in emergence of resistance. Based on the results of Monte Carlo simulation, the optimal dose of MFX to achieve an adequate antimicrobial efficacy should be 10 mg/kg, but it is important take into account that fluoroquinolones are substrates of efflux pumps, and this fact may determine that assumption of linear pharmacokinetics at high doses of MFX may be incorrect.     

金融市场随机波动：基于文献综述的视角

探讨了金融市场收益率存在历史、隐含和现实的三类随机波动现象,并呈现出尖峰厚尾、杠杆、集群、微笑、溢出、长记忆、信息流、共生波动等分布特征。进一步归纳了基于不同分布特征的随机波动的GARCH、SV、制度转换、阀值模型等模型,梳理出重点SV模型的三类估计方法：基于矩法、极大似然法估计和马尔科夫链蒙特卡罗（MCMC）方法,以及有效估计SV模型后,其在收益波动率预测、风险管理上的应用。