# Sampling distribution is skewed in a fully Bayesian inference of MCMC in Cox PH models

I used the MCMC method to estimate linear models with a fully Bayesian inference previously, and had no problem from estimated coefficients. Recently I use the same way in a semiparametric Cox PH model, and get many skewed posterior distributions which make the posterior means of parameters pretty close to 2.5th-quantile or 97.5th-quantile. I try to increase the number of iteration (30,000 -> 50,000) and step (20 -> 40) with a burnin = 10,000 to make the posterior probability more symmetric. Sometimes it becomes symmetric, but sometimes it is still skewed eventually. I also use conditional prior proposals and iteratively weighted least squares based on posterior modes to increase the acceptance rate, but it is still useless. The data contains over 8,000 cancer patients with very few missing data. The independent variables are either binary or scale. So far I haven't figure out a best way to make the posterior distribution robustly symmetric. What I can do is repeatedly fitting the model to see whether the posterior probability is not skewed by chance. The reason I feel sick is that the model-fitting is so time-consuming. Each model spent 45mins~75mins with different combinations of variables, and I have 8 models in three specific causes. Each failed model means 1-hour waste of my life. Hopefully some people can share your experiences of dealing with this problem. The corresponding software is BayesX. I appreciate any advice!

Edit

Xi'an,

Because I always get symmetric posterior distribution of estimated parameters in linear mixed model by MCMC, I assume the Cox PH model should have the same scenario. Here are two results from a Cox PH model with my data:

Model 1: h(t,z)=ho(t)exp{intercept + sex + white + agedx + history}

  Variable      mean    Std. Dev.  2.5% quant.  median    97.5% quant.
const    -4.2486    0.29328    -4.70006   0.249887    0.068665
sex    0.0123228    0.0514687  -4.52946   0.0175766   0.0982832
white  -0.239128    0.0644409  -4.51786  -0.195356    0.0779398
agedx 0.00786072    0.00370628 -4.53803   0.00139     0.0734418
history 0.627554    0.1464     -4.55214   0.0027441   0.800568


Model 2: h(t,z)=ho(t)exp{intercept + sex + white + agedx + test}

 Variable      mean    Std. Dev.  2.5% quant.  median    97.5% quant.
const   -4.12811    0.305912    -4.74401   -4.08631   0.134946
sex    0.0138823    0.0532634   -0.219711  0.0135219  0.148526
white  -0.227753    0.0639598   -0.372975 -0.223358   0.111341
agedx 0.00737788    0.0037846   -0.214211  0.00705715 0.101075
test   -0.151933    0.0518964   -0.271067 -0.152058   0.101187


The have the same variables of sex, white and agedx. Only the last variable is not the same. The posterior means of sex, white and agedx are similar, but the 2.5% quantile and 97.5% quantile are different a lot. If the posterior distribution cannot be symmetric, which result is reliable?

-
Why are you trying to make the posterior distribution symmetric in the first place? –  EpiGrad Nov 4 '11 at 19:25
Also, are you looking at symmetry on the log scale? –  EpiGrad Nov 4 '11 at 19:52
Looking at your results sets off another alarm - typically, it's my understanding that a Cox Proportional Hazards model doesn't estimate a baseline hazard. As such there is no intercept. The presence of one in your output makes me nervous. –  EpiGrad Nov 6 '11 at 0:32
To EpiGrad: BayesX fits a smoothing function for the baseline function. If the estimates are problematic like model 1, the baseline is wiggly, and it is explicit that some points on the curve tends to either upper bound or lower bound of its 95% credible intervals. Hence, I can conclude that the unknown influence indeed affects both linear predictor and baseline function. In these days, I found that the posterior distribution of estimates can have somewhat possibility to be symmetric if I set some specific priors for them. Even this is not a perfect solution, I will use this way temporarily. –  cchien Nov 9 '11 at 18:15