The below parametric Weibull survival AFT (survival R package) and its Bayesian alternative (Alvares 2021) under low informative priors should result in relatively similar results but they do not.

Any ideas on why?

library(rjags)
library(survival)

mc<- structure(list(V1 = c(6.3, 6.8, 4.7, 6.7, 6, 6.3, 6.9, 6.9, 5.6, 
                        7.1, 8, 6.9, 6.1, 5.8, 6.6, 8.1, 6, 6.6, 5.6, 6.1), 
                 V2 = structure(c(1L, 
                                  1L, 1L, 2L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 
                                  2L, 2L, 2L), levels = c("0", "1"), class = "factor"), 
                 V3 = c(0, 
                        0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 0, 1, 1, 0, 0, 1), 
                 X = c(1.65, 
                       1.79, 1.77, 2.74, 4.62, 0.86, 2.7, 3.47, 1.95, 4.3, 2.52, 1.95, 
                       2.78, 2.7, 4.25, 0.19, 3.76, 0.68, 1.48, 0.06), 
                 Time = c(38.863, 
                          24.2, 120.373, 16.5, 115.533, 8.613, 106.997, 105.71, 105.963, 
                          35.42, 107.613, 14.113, 105.82, 14.377, 99.253, 34.793, 96.657, 
                          94.347, 91.19, 8.613), 
                 Status = c(0, 0, 0, 1, 0, 1, 0, 0, 0, 
                            1, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1)), 
            row.names = c(1L, 2L, 4L, 
                          5L, 7L, 9L, 11L, 12L, 14L, 20L, 21L, 24L, 25L, 31L, 32L, 33L, 
                          37L, 40L, 41L, 42L), class = "data.frame")

BayesSurv<-function(form="X", data=mc, model="weibull",
                    prior_mu=0, prior_prec=1,var,
                    eventvariable="Status",
                    timevariable="Time",
                    niter=10000, nthin=10,nburn=1000){
  .formula <- reformulate(form)
  data$event<-data[,eventvariable]
  data$time<-data[,timevariable]
  X <- model.matrix(object=.formula, data = data)
  # Uncensored observations
  time_un <- data$time[data$event==1]
  cen_un <- data$time[data$event==1]
  is.censored_un <- rep(0, length(cen_un))
  X_un <- X[data$event==1, ]
  # Right-censored observations
  cen_cen <- data$time[data$event==0]
  time_cen <- rep(NA, length(cen_cen))
  is.censored_cen <- rep(1, length(cen_cen))
  X_cen <- X[data$event==0, ]
  # Uncensored and right-censored observations
  time <- c(time_un, time_cen)
  cen <- c(cen_un, cen_cen)
  is.censored <- c(is.censored_un, is.censored_cen)
  Xnew <- rbind(X_un, X_cen)
  
  cat("model{ #Uncensored and right−censored observations
for(i in 1:n){
is.censored[i]~dinterval(time[i],cen[i])
time[i]~dweib(alpha,lambda[i])
lambda[i]<-exp(-mu[i]*alpha)
mu[i]<-inprod(beta[],X[i,])}

#Prior distributions
for(l in 1:(Nbetas)){
beta[l]~dnorm(mu.betaMain[l],prec.betaMain[l])}
alpha ~ dunif(0, 10)
}",file = "AFT.txt")
  
  d.jags <- list(n = nrow(X), time = time, cen = cen, X = Xnew, is.censored = is.censored, Nbetas = ncol(X),mu.betaMain=prior_mu,prec.betaMain=prior_prec)
  i.jags <- function() list(beta = rnorm(ncol(X),sd = 1), alpha = runif(1,min = 0.1,max = 10))
  p.jags <- c("beta", "alpha")
  
  library("rjags")
  m1 <- rjags::jags.model(data = d.jags, file = "AFT.txt", inits = i.jags, n.chains = 1)
  
  # Burn-in and update
  update(m1, nburn)
  res <- coda.samples(m1, variable.names = p.jags, n.iter = niter, thin = nthin)
  results <- as.mcmc(do.call(rbind, res))
  
  return(results)
}


freqwei<-summary(survreg(Surv(Time,Status)~V1 + V2 + V3 + X, data=mc))


outn<-BayesSurv(form="V1 + V2 + V3 + X",prior_mu=c(0,rep(0,5)),prior_prec=c(rep(0.001,5)),
                           data = mc,niter = 20000)
bayewei=summary(outn)
bayewei

BOut<-cbind(freqwei$coefficients,
        bayewei[2]$quantiles[,3][-1],
        bayewei[1]$statistics[,1][-1])
BOut
    

………………………………………………………………………………………