Bayesian statistics in health economic evaluations
Bayes theoremHello folks, I hope you had some break time during summer as surely I did! After a whole year of stress and work it was nice to have some vacation period and to clear my mind for a while. Now that I am back to work I feel recharged and I am ready for a new year. One of my objective for this academic year is to find more time to dedicate to some new research projects as last year I only managed to do very little as most of my research energies went into the writing up of a research grant application. This year I hope to find more time to do something different, at least research wise.
So, with that spirit in mind, let’s start from today’s post where I follow-up from a past post introducing the concept of how to perform economic evaluations using standard statistical methods and power it up to what I normally do in this field, use Bayesian statistics! Perhaps some of you will not believe me but over time I am really sure that using Bayesian statistics made my life much easier when coming down to fit relatively complex models to health economics data. Since nowadays this seems to be very common in the literature, there is even more reason to go fully Bayesian when doing these analyses as the degree of flexibility it grants is so much more compared to what standard methods can typically achieve. Of course this is the opinion of someone totally biased! But before raising your finger, please try to come to the end of this post.
As I already mentioned in previous posts, the usual analysis task in economic evaluation based on individual-level data (e.g. QALYs and Total costs computed over a trial period) can be quite challenging due to the presence of a series of complexities that affect the data that need to be taken into account when choosing the statistical methods to use in the analysis; examples include: correlation between effects and costs, skewness of the outcome data, presence of structural values in the data. We saw before that different types of methods exist to deal with each of these problems but the general challenge comes from the fact that often these elements are present jointly in a single dataset and therefore the different methods used to handle each of them need to be combined in some way to perform the analysis. This, however, is easier said than done since, particularly under a frequentist framework, the complexity of fitting all these methods in combination with the need to quantify the impact of uncertainty on the results (e.g. via bootstrapping methods) can lead to extremely difficult-to-fit or expensive-to-implement models. This I believe the key reason why analysts often pretend to ignore some of these problems and prefer to implement easier methods in the hope that results will not be too much affected. Despite understanding their point, I feel that if they knew models that can account for all these problems together, then they would also like to fit them to improve the reliability of the results they obtain. Well, that is why today I talk about fitting the model under a Bayesian framework, which allows to achieve this task at the cost of learning a bit about Bayesian inference and how to interpret it.
Let’s start by simulating some non-Normal bivariate cost and QALY data from an hypothetical study for a total of \(300\) patients assigned to two competing intervention groups (\(t=0,1\)). When generating the data, we can try to mimic the typical skewness features of the outcome data by using alternative distributions such as Gamma for costs and Beta for QALYs. We also generate indicator variables that are used in order to determine which individuals should be assigned “structural values”, namely zero costs and one QALYs. The proportions of individuals assigned to these values is obtained by setting the probability of the Bernoulli distribution used to create the indicators.
set.seed(768)
n <- 300
id <- seq(1:n)
trt <- c(rep(0, n/2),rep(1, n/2))
mean_e1 <- c(0.5)
mean_e2 <- c(0.7)
sigma_e <- 0.15
tau1_e <- ((mean_e1*(1-mean_e1))/(sigma_e^2)-1)
tau2_e <- ((mean_e2*(1-mean_e2))/(sigma_e^2)-1)
alpha1_beta <- tau1_e*mean_e1
beta1_beta <- tau1_e*(1-mean_e1)
alpha2_beta <- tau2_e*mean_e2
beta2_beta <- tau2_e*(1-mean_e2)
e1 <- rbeta(n/2, alpha1_beta, beta1_beta)
e2 <- rbeta(n/2, alpha2_beta, beta2_beta)
mean_c1 <- 500
mean_c2 <- 1000
sigma_c <- 300
tau1_c <- mean_c1/(sigma_c^2)
tau2_c <- mean_c2/(sigma_c^2)
ln.mean_c1 <- log(500) + 5*(e1-mean(e1))
c1 <- rgamma(n/2, (exp(ln.mean_c1)/sigma_c)^2, exp(ln.mean_c1)/(sigma_c^2))
ln.mean_c2 <- log(1000) + 5*(e2-mean(e2)) + rgamma(n/2,0,tau2_c)
c2 <- rgamma(n/2, (exp(ln.mean_c2)/sigma_c)^2, exp(ln.mean_c2)/(sigma_c^2))
QALYs <- c(e1,e2)
Costs <- c(c1,c2)
p_zeros <- 0.25
d_zeros <- rbinom(n, 1, p_zeros)
p_ones <- 0.25
d_ones <- rbinom(n, 1, p_ones)
QALYs <- ifelse(d_ones==1,1,QALYs)
Costs <- ifelse(d_zeros==1,0,Costs)
data_sim_ec <- data.frame(id, trt, QALYs, Costs, d_zeros, d_ones)
data_sim_ec <- data_sim_ec[sample(1:nrow(data_sim_ec)), ]In the code above, after simulating QALY and Cost data using Beta and Gamma distribution, indicator variables for the zero and one values were generated for each individual in the sample from a Bernoulli distribution. Whenever the indicator takes value 1, it denotes the presence of a structural value and the corresponding outcome value is then set equal to zero (Costs) or one (QALYs). We can now compute the correlation between variables and plot the two outcome variables against each other to show how the presence of these structural values affect their corresponding association pattern.
#empirical correlation between e and c (across groups)
cor(data_sim_ec$QALYs,data_sim_ec$Costs)# [1] 0.3582116#scatterplot of e and c data by group
library(ggplot2)
data_sim_ec$trtf <- factor(data_sim_ec$trt)
levels(data_sim_ec$trtf) <- c("old","new")
ggplot(data_sim_ec, aes(x=QALYs, y=Costs)) +
geom_point(size=2, shape=16) + theme_classic() +
facet_wrap(~trtf)
In addition, we can also produce histograms of the distribution of the outcomes by treatment group to have a rough idea of the amount of structural values by type of outcome and treatment group in our sample.
data_sim_ec$trtf <- factor(data_sim_ec$trt)
levels(data_sim_ec$trtf) <- c("old","new")
QALY_hist <- ggplot(data_sim_ec, aes(x=QALYs))+
geom_histogram(color="black", fill="grey")+
facet_grid(trtf ~ .) + theme_classic()
Tcost_hist <- ggplot(data_sim_ec, aes(x=Costs))+
geom_histogram(color="black", fill="grey")+
facet_grid(trtf ~ .) + theme_classic()
gridExtra::grid.arrange(QALY_hist, Tcost_hist, nrow = 1, ncol = 2)
Step 1: fit a standard normal model
In order to explain the basics of how to fit a Bayesian model, let’s start by considering a (kind of) standard model based on Normal distributions for both QALYs and Total costs. However, I will slightly modify the model to allow for the correlation between the two outcomes, that is we fit a bivariate normal model \(p(e,c\mid \boldsymbol \theta)\), where \(e\) and \(c\) denote the QALYs and Total cost variables measured for each patient in the trial while \(\boldsymbol \theta\) denote the set of parameters indexing the model, including the key quantities of interest for the economic evaluations, i.e. the treatment-specific mean effect and cost \(\mu_{et}\) and \(\mu_{ct}\). To ease the task of modelling the data, we can re-express the joint distribution as:
\[p(e,c\mid \boldsymbol \theta) = p(e\mid \boldsymbol \theta_e) p(c \mid e, \boldsymbol \theta_c), \]
where \(p(e\mid \boldsymbol \theta_e)\) is the marginal distribution of the effects and \(p(c \mid e, \boldsymbol \theta_c)\) is the conditional distribution of the cost given the effects, each indexed by corresponding set of parameters. The main reason for factoring the joint distribution into this product is the possibility to specify univariate distributions for \(e\) and \(c\), rather than a single bivariate distribution. This can be helpful when, for example, different covariates are considered for the two outcomes as it allows a higher degree of flexibility in specifying the model for each variable. But how are we going to fit the model? well, for that we can rely on freely-available Bayesian software which allows model fitting in a relatively simple way at the cost of learning how to code up the model in this new language. In this post I will consider the JAGS software although this is only one of the many that can be used. For the sake of making things clearer I will not focus here on the details of how the software works and which types of algorithms it uses to implement the model, but I will jump straight into the coding part. First, we need to write the code of the model into a txt file that will then be read by the program after providing the data as input. We can do all this in R.
model_bn <- "
model {
#model specification
for(i in 1:n){
QALYs[i] ~ dnorm(nu_e[i],tau_e)
nu_e[i] <- beta0 + beta1*trt[i]
Costs[i] ~ dnorm(nu_c[i],tau_c)
nu_c[i] <- gamma0 + gamma1*trt[i] + gamma2*QALYs[i]
}
#prior specification
tau_e <- 1/ss_e
ss_e <- s_e*s_e
tau_c <- 1/ss_c
ss_c <- s_c*s_c
s_c ~ dunif(0,1000)
s_e ~ dunif(0,1000)
beta0 ~ dnorm(0,0.000001)
beta1 ~ dnorm(0,0.000001)
gamma0 ~ dnorm(0,0.000001)
gamma1 ~ dnorm(0,0.000001)
gamma2 ~ dnorm(0,0.000001)
}
"
writeLines(model_bn, con = "model_bn.txt")In the code above I first specify the model structure, i.e. assign normal distributions to QLAYs and Costs variable indexed by two parameters, the means \(\nu\) and precisions \(\tau\) (note that precisions correspond to inverse of the variance \(\tau=1/\sigma^2\)) since these are the default parameters used by JAGS to specify a normal distribution. For each outcome then I specify the mean structure, i.e. the mean of \(e\) depennds only on the treatment indicator while the mean of \(c\) depend both on treatment indicator and \(e\) (this is a conditional cost model!). Next, I specify the priors for non-deterministic parameters, namely using uniform distributions for standard deviations and normal distributions for regression coefficients. Finally, I save the model as a txt file in the current wd using the writeLines function. The model is now written and we can fit it by calling the JAGS software directly from R through dedicated functions. Before that, we need to convert the data as input for the software. Then, we load the package R2jags which allows to call the software from R through the function jags and after providing some technical parameters needed to run the model.
#save data input
n <- dim(data_sim_ec)[1]
QALYs <- data_sim_ec$QALYs
Costs <- data_sim_ec$Costs
trt <- data_sim_ec$trt
#load package and provide algorithm parameters
library(R2jags)
set.seed(2345) #set seed for reproducibility
datalist<-list("n","QALYs","Costs","trt") #pass data into a list
#set up initial values for algorithm
inits1 <- list(.RNG.name = "base::Wichmann-Hill", .RNG.seed = 1)
inits2 <- list(.RNG.name = "base::Wichmann-Hill", .RNG.seed = 2)
#set parameter easimates to save
params<-c("beta0","beta1","gamma0","gamma1","gamma2","s_c","s_e","nu_c","nu_e")
filein<-"model_bn.txt" #name of model file
n.iter<-20000 #n of iterationsWe are now ready to fit the model, whic we can do by typing
#fit model
jmodel_bn<-jags(data=datalist,inits=list(inits1,inits2),
parameters.to.save=params,model.file=filein,
n.chains=2,n.iter=n.iter,n.thin=1)# Compiling model graph
# Resolving undeclared variables
# Allocating nodes
# Graph information:
# Observed stochastic nodes: 600
# Unobserved stochastic nodes: 7
# Total graph size: 1345
#
# Initializing modelAfter some time needed for the model to run, we end up with something like this
#posterior results
print(jmodel_bn)# Inference for Bugs model at "model_bn.txt", fit using jags,
# 2 chains, each with 20000 iterations (first 10000 discarded)
# n.sims = 20000 iterations saved
# mu.vect sd.vect 2.5% 25% 50% 75% 97.5% Rhat
# beta0 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# beta1 0.116 0.027 0.063 0.098 0.116 0.134 0.168 1.001
# gamma0 -147.064 136.100 -414.191 -238.684 -146.655 -56.422 120.626 1.001
# gamma1 443.594 89.617 268.643 383.247 443.401 504.381 618.307 1.001
# gamma2 1003.566 185.780 642.338 878.610 1002.303 1128.245 1368.741 1.001
# nu_c[1] 653.790 99.155 458.610 587.817 653.583 719.540 847.442 1.001
# nu_c[2] 1019.968 62.799 897.609 977.709 1019.448 1062.243 1143.497 1.001
# nu_c[3] 345.232 68.708 211.853 298.661 344.974 391.561 480.094 1.001
# nu_c[4] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[5] 625.502 65.446 495.707 581.253 625.514 669.387 754.077 1.001
# nu_c[6] 1194.620 66.124 1065.903 1150.333 1194.031 1238.974 1323.582 1.001
# nu_c[7] 460.143 62.371 338.056 418.112 460.125 502.151 583.002 1.001
# nu_c[8] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[9] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[10] 563.317 62.574 440.263 521.301 563.189 605.310 686.902 1.001
# nu_c[11] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[12] 1191.867 65.950 1063.432 1147.673 1191.312 1235.984 1320.751 1.001
# nu_c[13] 994.028 63.705 870.006 951.189 993.596 1036.803 1120.059 1.001
# nu_c[14] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[15] 1165.410 64.460 1039.591 1122.082 1164.928 1208.722 1291.399 1.001
# nu_c[16] 809.178 78.808 654.435 756.291 809.094 861.036 964.868 1.001
# nu_c[17] 340.903 69.063 206.924 294.139 340.636 387.481 476.537 1.001
# nu_c[18] 872.921 72.138 732.591 824.388 872.264 920.863 1015.917 1.001
# nu_c[19] 500.837 61.750 379.933 459.165 500.833 542.352 622.304 1.001
# nu_c[20] 973.213 64.682 847.399 929.984 972.492 1016.510 1100.975 1.001
# nu_c[21] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[22] 974.684 64.606 848.977 931.523 973.961 1017.918 1102.260 1.001
# nu_c[23] 837.468 75.693 690.179 786.761 837.409 887.411 987.207 1.001
# nu_c[24] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[25] 806.961 79.061 651.919 753.944 806.819 859.092 963.011 1.001
# nu_c[26] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[27] 361.790 67.419 230.804 316.064 361.775 407.072 494.146 1.001
# nu_c[28] 1024.197 62.685 901.741 981.990 1023.696 1066.397 1147.427 1.001
# nu_c[29] 383.746 65.890 255.443 338.978 383.877 428.377 513.959 1.001
# nu_c[30] 227.424 80.696 70.418 173.207 227.225 282.027 385.480 1.001
# nu_c[31] 929.086 67.436 797.581 883.929 928.507 974.166 1062.601 1.001
# nu_c[32] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[33] 823.970 77.151 673.371 772.360 823.975 874.825 976.251 1.001
# nu_c[34] 1132.331 63.085 1009.282 1089.397 1132.140 1174.378 1255.022 1.001
# nu_c[35] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[36] 1118.608 62.681 996.209 1076.318 1118.460 1160.549 1240.764 1.001
# nu_c[37] 359.446 67.594 228.237 313.609 359.371 404.935 492.136 1.001
# nu_c[38] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[39] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[40] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[41] 588.930 108.778 374.370 516.785 588.994 661.198 802.627 1.001
# nu_c[42] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[43] 455.381 62.503 332.926 413.332 455.433 497.558 578.525 1.001
# nu_c[44] 740.735 87.187 570.340 682.575 740.532 798.336 912.542 1.001
# nu_c[45] 430.231 63.398 306.488 387.377 430.247 473.506 555.669 1.001
# nu_c[46] 441.056 62.972 317.871 398.451 441.096 483.796 565.390 1.001
# nu_c[47] 1169.891 64.689 1043.692 1126.451 1169.458 1213.259 1296.198 1.001
# nu_c[48] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[49] 288.014 73.972 144.226 238.023 287.678 338.282 432.055 1.001
# nu_c[50] 928.271 67.496 796.571 883.110 927.761 973.420 1061.858 1.001
# nu_c[51] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[52] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[53] 1072.145 62.071 950.668 1030.223 1071.914 1113.477 1193.582 1.001
# nu_c[54] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[55] 988.362 63.950 863.777 945.424 987.777 1031.150 1114.910 1.001
# nu_c[56] 128.235 93.581 -53.831 65.024 127.986 190.842 312.124 1.001
# nu_c[57] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[58] 146.811 91.025 -29.725 85.618 146.778 208.016 325.036 1.001
# nu_c[59] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[60] 607.099 106.033 398.304 536.753 606.724 677.340 815.537 1.001
# nu_c[61] 152.831 90.209 -21.966 92.077 152.643 213.540 329.462 1.001
# nu_c[62] 501.200 61.748 380.299 459.536 501.200 542.677 622.629 1.001
# nu_c[63] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[64] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[65] 972.467 64.721 846.548 929.167 971.758 1015.811 1100.624 1.001
# nu_c[66] 556.390 62.378 433.787 514.490 556.439 598.405 679.073 1.001
# nu_c[67] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[68] 611.944 105.307 404.612 542.031 611.649 681.797 818.630 1.001
# nu_c[69] 1135.728 63.200 1012.346 1092.687 1135.406 1177.868 1258.411 1.001
# nu_c[70] 1055.716 62.140 934.313 1013.709 1055.420 1097.361 1177.248 1.001
# nu_c[71] 831.932 76.284 683.435 780.881 831.983 882.173 982.980 1.001
# nu_c[72] 778.404 82.440 617.205 723.291 778.203 832.481 940.706 1.001
# nu_c[73] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[74] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[75] 865.762 72.822 723.928 816.706 865.138 914.183 1010.231 1.001
# nu_c[76] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[77] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[78] 1084.918 62.121 963.388 1043.058 1084.629 1126.439 1206.235 1.001
# nu_c[79] 200.283 84.019 37.557 143.579 200.108 257.303 364.512 1.001
# nu_c[80] 331.705 69.843 196.626 284.343 331.304 378.963 469.184 1.001
# nu_c[81] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[82] 300.374 72.737 159.373 251.200 300.150 349.635 442.519 1.001
# nu_c[83] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[84] 1180.882 65.291 1053.657 1137.075 1180.348 1224.545 1308.583 1.001
# nu_c[85] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[86] 1035.128 62.435 912.887 993.050 1034.908 1077.069 1157.367 1.001
# nu_c[87] 1049.819 62.200 928.122 1007.790 1049.655 1091.573 1171.365 1.001
# nu_c[88] 480.979 61.940 359.345 439.130 480.732 522.537 602.747 1.001
# nu_c[89] 1129.820 63.003 1006.979 1087.061 1129.528 1171.899 1252.262 1.001
# nu_c[90] 972.579 64.715 846.671 929.290 971.854 1015.907 1100.674 1.001
# nu_c[91] 918.919 68.194 785.895 873.028 918.408 964.483 1053.651 1.001
# nu_c[92] 431.336 63.351 307.756 388.563 431.312 474.563 556.622 1.001
# nu_c[93] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[94] 449.700 62.676 326.853 407.428 449.703 492.178 573.204 1.001
# nu_c[95] 1060.214 62.106 938.643 1018.300 1059.854 1101.722 1181.800 1.001
# nu_c[96] 453.978 62.544 331.473 411.816 453.999 496.181 577.207 1.001
# nu_c[97] 775.863 82.750 614.126 720.477 775.750 830.194 939.044 1.001
# nu_c[98] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[99] 1075.338 62.075 953.857 1033.428 1075.046 1116.783 1197.246 1.001
# nu_c[100] 1265.739 71.736 1126.516 1217.434 1265.248 1314.108 1407.430 1.001
# nu_c[101] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[102] 567.378 62.701 444.136 525.268 567.360 609.454 691.153 1.001
# nu_c[103] 471.741 62.102 350.168 430.014 471.659 513.439 593.879 1.001
# nu_c[104] 390.904 65.439 263.424 346.660 391.057 435.257 520.019 1.001
# nu_c[105] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[106] 545.692 62.125 423.762 503.696 545.709 587.546 667.842 1.001
# nu_c[107] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[108] 1246.361 70.006 1110.454 1199.244 1245.694 1293.687 1384.062 1.001
# nu_c[109] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[110] 453.798 62.550 331.319 411.643 453.823 496.027 577.068 1.001
# nu_c[111] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[112] 829.187 76.581 679.831 777.979 829.247 879.665 980.609 1.001
# nu_c[113] 1173.521 64.882 1047.074 1129.969 1173.028 1216.966 1300.152 1.001
# nu_c[114] 579.849 63.146 455.542 537.246 579.870 622.206 704.546 1.001
# nu_c[115] 878.268 71.639 739.032 830.161 877.622 926.056 1020.183 1.001
# nu_c[116] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[117] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[118] 1087.045 62.138 965.397 1045.166 1086.754 1128.640 1208.452 1.001
# nu_c[119] 269.480 75.916 121.941 218.248 269.373 320.905 417.606 1.001
# nu_c[120] 217.486 81.893 58.804 162.411 217.173 272.976 377.719 1.001
# nu_c[121] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[122] 411.898 64.255 286.333 368.495 411.972 455.562 538.746 1.001
# nu_c[123] 340.754 69.076 206.779 293.997 340.491 387.344 476.440 1.001
# nu_c[124] 393.508 65.281 266.302 349.330 393.631 437.807 522.578 1.001
# nu_c[125] 993.860 63.712 869.841 951.039 993.414 1036.649 1119.920 1.001
# nu_c[126] 918.847 68.200 785.814 872.950 918.333 964.407 1053.593 1.001
# nu_c[127] 942.708 66.491 813.072 898.199 941.958 987.311 1074.899 1.001
# nu_c[128] 600.639 64.066 474.121 557.393 600.442 643.701 726.447 1.001
# nu_c[129] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[130] 238.751 79.362 84.821 185.310 238.637 292.281 394.532 1.001
# nu_c[131] 291.722 73.596 148.600 241.944 291.345 341.674 435.143 1.001
# nu_c[132] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[133] 272.240 75.620 125.369 221.126 272.136 323.435 419.896 1.001
# nu_c[134] 862.012 73.188 719.604 812.768 861.429 910.596 1007.066 1.001
# nu_c[135] 527.533 61.839 406.251 485.702 527.602 569.191 649.584 1.001
# nu_c[136] 1109.593 62.471 987.488 1067.352 1109.397 1151.373 1231.314 1.001
# nu_c[137] 357.166 67.767 225.565 311.261 357.093 402.750 490.291 1.001
# nu_c[138] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[139] 1014.176 62.971 891.107 971.770 1013.619 1056.509 1138.054 1.001
# nu_c[140] 1096.916 62.249 975.076 1054.895 1096.604 1138.436 1218.545 1.001
# nu_c[141] 837.029 75.740 689.575 786.289 836.974 886.966 986.871 1.001
# nu_c[142] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[143] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[144] 440.123 63.006 316.866 397.468 440.150 482.925 564.565 1.001
# nu_c[145] 401.262 64.828 274.856 357.348 401.464 445.279 529.296 1.001
# nu_c[146] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[147] 364.962 67.185 234.535 319.386 365.043 410.060 497.141 1.001
# nu_c[148] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[149] 463.035 62.297 341.267 420.976 463.007 504.950 585.569 1.001
# nu_c[150] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[151] 165.754 88.480 -5.480 106.133 165.719 225.351 339.038 1.001
# nu_c[152] 377.056 66.333 247.901 331.963 377.180 421.751 507.572 1.001
# nu_c[153] 507.359 61.735 386.562 465.664 507.384 548.947 628.816 1.001
# nu_c[154] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[155] 1029.794 62.549 907.183 987.624 1029.433 1071.778 1152.676 1.001
# nu_c[156] 676.339 95.935 487.461 612.554 676.129 739.855 863.658 1.001
# nu_c[157] 452.035 62.603 329.448 409.841 452.074 494.344 575.421 1.001
# nu_c[158] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[159] 1037.890 62.382 915.811 995.760 1037.673 1079.745 1159.808 1.001
# nu_c[160] 1097.436 62.257 975.554 1055.434 1097.111 1138.986 1219.053 1.001
# nu_c[161] 992.475 63.771 868.361 949.591 991.984 1035.225 1118.779 1.001
# nu_c[162] 197.107 84.419 33.531 140.025 197.024 254.196 361.972 1.001
# nu_c[163] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[164] 432.218 63.315 308.621 389.404 432.183 475.332 557.546 1.001
# nu_c[165] 604.453 64.258 477.469 561.030 604.294 647.567 730.810 1.001
# nu_c[166] 1108.516 62.448 986.372 1066.250 1108.262 1150.275 1230.153 1.001
# nu_c[167] 490.674 61.820 369.481 449.052 490.515 532.279 612.184 1.001
# nu_c[168] 297.238 73.045 155.364 247.821 296.911 346.776 439.597 1.001
# nu_c[169] 757.527 85.033 591.281 700.792 757.319 813.510 924.869 1.001
# nu_c[170] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[171] 769.315 83.557 606.126 713.484 769.219 824.229 933.923 1.001
# nu_c[172] 503.620 61.741 382.632 461.915 503.588 545.131 624.894 1.001
# nu_c[173] 130.676 93.242 -50.660 67.756 130.524 193.125 313.898 1.001
# nu_c[174] 959.564 65.437 832.146 915.742 958.767 1003.510 1089.361 1.001
# nu_c[175] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[176] 170.081 87.908 -0.093 110.791 169.977 229.324 342.003 1.001
# nu_c[177] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[178] 305.910 72.200 166.109 257.031 305.625 354.740 447.209 1.001
# nu_c[179] 481.343 61.934 359.708 439.493 481.098 522.917 603.097 1.001
# nu_c[180] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[181] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[182] 534.462 61.927 412.766 492.529 534.620 576.220 656.674 1.001
# nu_c[183] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[184] 1036.865 62.401 914.755 994.690 1036.621 1078.762 1158.836 1.001
# nu_c[185] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[186] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[187] 228.208 80.602 71.421 174.065 227.969 282.654 386.077 1.001
# nu_c[188] 453.709 62.552 331.241 411.562 453.718 495.945 576.989 1.001
# nu_c[189] 1079.333 62.088 958.114 1037.493 1078.952 1120.800 1200.873 1.001
# nu_c[190] 1072.557 62.071 951.115 1030.611 1072.321 1113.916 1193.961 1.001
# nu_c[191] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[192] 436.323 63.151 313.002 393.679 436.323 479.210 561.095 1.001
# nu_c[193] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[194] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[195] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[196] 1017.705 62.864 895.085 975.433 1017.095 1060.072 1141.432 1.001
# nu_c[197] 162.162 88.957 -10.108 102.143 162.094 222.163 335.939 1.001
# nu_c[198] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[199] 1067.667 62.075 945.754 1025.813 1067.496 1109.081 1188.903 1.001
# nu_c[200] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[201] 323.111 70.601 186.414 275.174 322.813 371.022 462.002 1.001
# nu_c[202] 1100.266 62.299 978.213 1058.213 1100.001 1141.787 1221.926 1.001
# nu_c[203] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[204] 923.430 67.853 790.929 877.931 922.885 968.765 1057.648 1.001
# nu_c[205] 808.787 78.852 653.982 755.927 808.727 860.683 964.561 1.001
# nu_c[206] 73.696 101.388 -123.777 5.095 73.187 141.186 273.770 1.001
# nu_c[207] 148.163 90.841 -27.996 87.123 148.121 209.274 325.924 1.001
# nu_c[208] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[209] 780.940 82.132 620.262 726.041 780.662 834.767 942.495 1.001
# nu_c[210] 772.066 83.217 609.404 716.391 771.970 826.645 936.222 1.001
# nu_c[211] 566.636 62.677 443.372 524.561 566.575 608.673 690.331 1.001
# nu_c[212] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[213] 1014.910 62.948 891.849 972.478 1014.298 1057.259 1138.751 1.001
# nu_c[214] 712.363 90.953 533.536 651.553 712.275 772.596 890.455 1.001
# nu_c[215] 992.563 63.767 868.445 949.670 992.060 1035.320 1118.845 1.001
# nu_c[216] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[217] 259.008 77.060 109.562 207.075 258.796 310.995 409.677 1.001
# nu_c[218] 910.403 68.862 776.619 864.265 909.869 956.347 1046.502 1.001
# nu_c[219] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[220] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[221] 1135.591 63.195 1012.223 1092.577 1135.274 1177.724 1258.271 1.001
# nu_c[222] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[223] 314.743 71.366 176.571 266.342 314.410 363.168 454.796 1.001
# nu_c[224] 166.929 88.324 -4.136 107.394 166.891 226.387 339.874 1.001
# nu_c[225] 1137.322 63.256 1013.909 1094.325 1137.065 1179.437 1260.155 1.001
# nu_c[226] 902.146 69.538 766.438 855.432 901.529 948.629 1039.498 1.001
# nu_c[227] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[228] 339.275 69.199 205.103 292.389 339.044 385.902 475.157 1.001
# nu_c[229] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[230] 1183.029 65.416 1055.654 1139.089 1182.527 1226.752 1310.951 1.001
# nu_c[231] 1048.167 62.221 926.479 1006.160 1047.966 1089.878 1169.762 1.001
# nu_c[232] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[233] 471.388 62.109 349.843 429.703 471.298 513.122 593.505 1.001
# nu_c[234] 756.718 85.136 590.368 699.890 756.510 812.821 924.330 1.001
# nu_c[235] 1015.861 62.919 892.961 973.494 1015.225 1058.222 1139.715 1.001
# nu_c[236] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[237] 433.514 63.262 310.113 390.746 433.533 476.546 558.628 1.001
# nu_c[238] 354.821 67.947 222.710 308.802 354.675 400.591 488.379 1.001
# nu_c[239] 1191.369 65.919 1062.956 1147.190 1190.825 1235.472 1320.154 1.001
# nu_c[240] 1008.980 63.140 885.696 966.597 1008.513 1051.395 1133.476 1.001
# nu_c[241] 710.029 91.269 530.568 649.105 709.985 770.399 888.743 1.001
# nu_c[242] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[243] 379.830 66.147 251.073 334.858 380.032 424.445 510.398 1.001
# nu_c[244] 472.134 62.094 350.590 430.446 472.051 513.847 594.245 1.001
# nu_c[245] 568.522 62.739 445.109 526.389 568.511 610.636 692.440 1.001
# nu_c[246] 525.065 118.679 290.495 445.865 525.356 604.315 757.512 1.001
# nu_c[247] 355.505 67.895 223.552 309.538 355.419 401.219 488.949 1.001
# nu_c[248] 1224.452 68.225 1091.572 1178.466 1223.865 1270.420 1358.569 1.001
# nu_c[249] 307.717 72.027 168.248 259.022 307.392 356.456 448.801 1.001
# nu_c[250] 1093.948 62.210 972.201 1052.079 1093.639 1135.454 1215.576 1.001
# nu_c[251] 313.416 71.489 174.968 264.926 313.081 361.911 453.697 1.001
# nu_c[252] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[253] 217.555 81.884 58.885 162.485 217.240 273.041 377.784 1.001
# nu_c[254] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[255] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[256] 274.688 75.359 128.526 223.717 274.514 325.774 421.700 1.001
# nu_c[257] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[258] 577.648 63.062 453.270 535.044 577.676 619.974 702.294 1.001
# nu_c[259] 288.377 73.935 144.676 238.407 288.093 338.608 432.318 1.001
# nu_c[260] 962.194 65.285 834.977 918.512 961.499 1005.948 1091.694 1.001
# nu_c[261] 289.305 73.841 145.714 239.348 289.004 339.481 433.145 1.001
# nu_c[262] 1081.546 62.099 960.130 1039.688 1081.260 1123.012 1202.804 1.001
# nu_c[263] 188.313 85.538 22.210 130.477 188.244 246.045 355.528 1.001
# nu_c[264] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[265] 56.582 103.919 -145.814 -14.036 56.083 125.730 261.937 1.001
# nu_c[266] 793.987 80.571 636.301 739.769 793.778 847.068 952.291 1.001
# nu_c[267] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[268] 1144.670 63.532 1020.502 1101.573 1144.416 1187.147 1268.349 1.001
# nu_c[269] 231.478 80.215 75.627 177.611 231.251 285.529 388.696 1.001
# nu_c[270] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[271] 327.033 70.252 191.121 279.389 326.818 374.628 465.338 1.001
# nu_c[272] 946.117 66.268 817.019 901.739 945.340 990.589 1077.745 1.001
# nu_c[273] 573.242 62.900 449.396 530.871 573.247 615.448 697.549 1.001
# nu_c[274] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[275] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[276] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[277] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[278] 1232.160 68.830 1098.163 1185.729 1231.552 1278.396 1367.565 1.001
# nu_c[279] 1009.228 63.132 885.931 966.839 1008.768 1051.619 1133.730 1.001
# nu_c[280] 164.753 88.613 -6.774 105.008 164.741 224.484 338.222 1.001
# nu_c[281] 1300.096 75.126 1155.050 1249.149 1299.742 1350.649 1449.069 1.001
# nu_c[282] 626.000 65.477 496.128 581.724 626.001 669.931 754.663 1.001
# nu_c[283] 257.657 77.210 107.917 205.546 257.496 309.786 408.557 1.001
# nu_c[284] 372.097 66.675 242.342 326.834 372.134 416.931 503.345 1.001
# nu_c[285] 715.856 90.481 537.992 655.407 715.729 775.835 893.483 1.001
# nu_c[286] 403.727 64.691 277.568 359.971 403.867 447.772 531.564 1.001
# nu_c[287] 1152.203 63.844 1027.452 1109.080 1151.926 1195.022 1276.850 1.001
# nu_c[288] 1004.235 63.307 880.782 961.825 1003.827 1046.725 1129.082 1.001
# nu_c[289] 426.669 63.551 302.407 383.738 426.792 469.977 552.243 1.001
# nu_c[290] 782.928 81.891 622.579 728.116 782.675 836.602 944.040 1.001
# nu_c[291] 95.061 98.280 -96.099 28.409 94.726 160.560 288.836 1.001
# nu_c[292] 398.652 64.977 271.841 354.557 398.829 442.788 526.913 1.001
# nu_c[293] 1142.963 63.466 1019.007 1099.823 1142.655 1185.382 1266.339 1.001
# nu_c[294] 397.086 65.069 270.034 352.930 397.223 441.330 525.547 1.001
# nu_c[295] 416.723 64.013 291.533 373.423 416.783 460.279 543.401 1.001
# nu_c[296] 155.756 89.815 -18.145 95.239 155.576 216.315 331.500 1.001
# nu_c[297] 988.642 63.937 864.104 945.693 988.055 1031.457 1115.224 1.001
# nu_c[298] 139.710 91.995 -38.810 77.675 139.478 201.461 320.181 1.001
# nu_c[299] 856.502 89.274 679.911 796.575 856.125 916.588 1031.153 1.001
# nu_c[300] 1116.904 62.638 994.629 1074.631 1116.739 1158.830 1238.954 1.001
# nu_e[1] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[2] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[3] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[4] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[5] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[6] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[7] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[8] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[9] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[10] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[11] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[12] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[13] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[14] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[15] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[16] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[17] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[18] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[19] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[20] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[21] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[22] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[23] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[24] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[25] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[26] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[27] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[28] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[29] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[30] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[31] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[32] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[33] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[34] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[35] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[36] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[37] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[38] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[39] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[40] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[41] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[42] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[43] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[44] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[45] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[46] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[47] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[48] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[49] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[50] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[51] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[52] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[53] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[54] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[55] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[56] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[57] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[58] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[59] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[60] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[61] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[62] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[63] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[64] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[65] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[66] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[67] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[68] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[69] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[70] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[71] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[72] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[73] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[74] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[75] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[76] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[77] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[78] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[79] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[80] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[81] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[82] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[83] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[84] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[85] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[86] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[87] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[88] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[89] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[90] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[91] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[92] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[93] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[94] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[95] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[96] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[97] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[98] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[99] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[100] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[101] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[102] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[103] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[104] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[105] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[106] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[107] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[108] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[109] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[110] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[111] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[112] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[113] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[114] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[115] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[116] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[117] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[118] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[119] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[120] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[121] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[122] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[123] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[124] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[125] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[126] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[127] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[128] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[129] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[130] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[131] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[132] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[133] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[134] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[135] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[136] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[137] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[138] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[139] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[140] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[141] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[142] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[143] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[144] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[145] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[146] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[147] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[148] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[149] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[150] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[151] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[152] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[153] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[154] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[155] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[156] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[157] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[158] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[159] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[160] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[161] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[162] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[163] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[164] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[165] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[166] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[167] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[168] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[169] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[170] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[171] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[172] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[173] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[174] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[175] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[176] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[177] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[178] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[179] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[180] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[181] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[182] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[183] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[184] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[185] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[186] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[187] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[188] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[189] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[190] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[191] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[192] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[193] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[194] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[195] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[196] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[197] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[198] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[199] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[200] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[201] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[202] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[203] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[204] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[205] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[206] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[207] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[208] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[209] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[210] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[211] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[212] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[213] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[214] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[215] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[216] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[217] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[218] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[219] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[220] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[221] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[222] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[223] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[224] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[225] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[226] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[227] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[228] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[229] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[230] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[231] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[232] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[233] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[234] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[235] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[236] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[237] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[238] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[239] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[240] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[241] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[242] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[243] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[244] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[245] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[246] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[247] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[248] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[249] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[250] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[251] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[252] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[253] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[254] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[255] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[256] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[257] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[258] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[259] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[260] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[261] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[262] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[263] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[264] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[265] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[266] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[267] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[268] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[269] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[270] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[271] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[272] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[273] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[274] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[275] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[276] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[277] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[278] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[279] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[280] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[281] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[282] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[283] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[284] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[285] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[286] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[287] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[288] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[289] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[290] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[291] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[292] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[293] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[294] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[295] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[296] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[297] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# nu_e[298] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[299] 0.656 0.019 0.619 0.643 0.656 0.669 0.693 1.001
# nu_e[300] 0.772 0.019 0.735 0.759 0.772 0.785 0.809 1.001
# s_c 759.603 31.405 701.198 737.789 758.687 780.114 824.233 1.001
# s_e 0.231 0.009 0.214 0.225 0.231 0.237 0.251 1.001
# deviance 4800.101 3.774 4794.804 4797.326 4799.423 4802.136 4809.117 1.001
# n.eff
# beta0 20000
# beta1 20000
# gamma0 20000
# gamma1 20000
# gamma2 20000
# nu_c[1] 13000
# nu_c[2] 6900
# nu_c[3] 9900
# nu_c[4] 12000
# nu_c[5] 12000
# nu_c[6] 8800
# nu_c[7] 9500
# nu_c[8] 20000
# nu_c[9] 12000
# nu_c[10] 11000
# nu_c[11] 12000
# nu_c[12] 8800
# nu_c[13] 7000
# nu_c[14] 20000
# nu_c[15] 8200
# nu_c[16] 9300
# nu_c[17] 9900
# nu_c[18] 8100
# nu_c[19] 9700
# nu_c[20] 7100
# nu_c[21] 12000
# nu_c[22] 7100
# nu_c[23] 8700
# nu_c[24] 20000
# nu_c[25] 9300
# nu_c[26] 20000
# nu_c[27] 9800
# nu_c[28] 6900
# nu_c[29] 9700
# nu_c[30] 12000
# nu_c[31] 7400
# nu_c[32] 20000
# nu_c[33] 9000
# nu_c[34] 7600
# nu_c[35] 20000
# nu_c[36] 7500
# nu_c[37] 9800
# nu_c[38] 12000
# nu_c[39] 20000
# nu_c[40] 12000
# nu_c[41] 14000
# nu_c[42] 20000
# nu_c[43] 9400
# nu_c[44] 11000
# nu_c[45] 9500
# nu_c[46] 9400
# nu_c[47] 8300
# nu_c[48] 12000
# nu_c[49] 11000
# nu_c[50] 7400
# nu_c[51] 12000
# nu_c[52] 20000
# nu_c[53] 7000
# nu_c[54] 20000
# nu_c[55] 7000
# nu_c[56] 15000
# nu_c[57] 20000
# nu_c[58] 14000
# nu_c[59] 20000
# nu_c[60] 14000
# nu_c[61] 14000
# nu_c[62] 9700
# nu_c[63] 20000
# nu_c[64] 20000
# nu_c[65] 7100
# nu_c[66] 10000
# nu_c[67] 20000
# nu_c[68] 14000
# nu_c[69] 7700
# nu_c[70] 7000
# nu_c[71] 8800
# nu_c[72] 9900
# nu_c[73] 20000
# nu_c[74] 20000
# nu_c[75] 8300
# nu_c[76] 20000
# nu_c[77] 20000
# nu_c[78] 7100
# nu_c[79] 13000
# nu_c[80] 9900
# nu_c[81] 12000
# nu_c[82] 10000
# nu_c[83] 12000
# nu_c[84] 8500
# nu_c[85] 12000
# nu_c[86] 6900
# nu_c[87] 6900
# nu_c[88] 9500
# nu_c[89] 7600
# nu_c[90] 7100
# nu_c[91] 7500
# nu_c[92] 9500
# nu_c[93] 12000
# nu_c[94] 9400
# nu_c[95] 7000
# nu_c[96] 9400
# nu_c[97] 9900
# nu_c[98] 20000
# nu_c[99] 7100
# nu_c[100] 11000
# nu_c[101] 12000
# nu_c[102] 11000
# nu_c[103] 9500
# nu_c[104] 9600
# nu_c[105] 12000
# nu_c[106] 10000
# nu_c[107] 12000
# nu_c[108] 10000
# nu_c[109] 20000
# nu_c[110] 9400
# nu_c[111] 20000
# nu_c[112] 8900
# nu_c[113] 8400
# nu_c[114] 11000
# nu_c[115] 8100
# nu_c[116] 20000
# nu_c[117] 12000
# nu_c[118] 7100
# nu_c[119] 11000
# nu_c[120] 12000
# nu_c[121] 12000
# nu_c[122] 9500
# nu_c[123] 9900
# nu_c[124] 9600
# nu_c[125] 7000
# nu_c[126] 7500
# nu_c[127] 7300
# nu_c[128] 12000
# nu_c[129] 20000
# nu_c[130] 12000
# nu_c[131] 11000
# nu_c[132] 20000
# nu_c[133] 11000
# nu_c[134] 8300
# nu_c[135] 10000
# nu_c[136] 7400
# nu_c[137] 9800
# nu_c[138] 12000
# nu_c[139] 6900
# nu_c[140] 7200
# nu_c[141] 8700
# nu_c[142] 12000
# nu_c[143] 20000
# nu_c[144] 9400
# nu_c[145] 9600
# nu_c[146] 12000
# nu_c[147] 9800
# nu_c[148] 12000
# nu_c[149] 9500
# nu_c[150] 12000
# nu_c[151] 13000
# nu_c[152] 9700
# nu_c[153] 9700
# nu_c[154] 20000
# nu_c[155] 6900
# nu_c[156] 12000
# nu_c[157] 9400
# nu_c[158] 12000
# nu_c[159] 6900
# nu_c[160] 7200
# nu_c[161] 7000
# nu_c[162] 13000
# nu_c[163] 12000
# nu_c[164] 9400
# nu_c[165] 12000
# nu_c[166] 7300
# nu_c[167] 9600
# nu_c[168] 10000
# nu_c[169] 10000
# nu_c[170] 12000
# nu_c[171] 10000
# nu_c[172] 9700
# nu_c[173] 14000
# nu_c[174] 7200
# nu_c[175] 12000
# nu_c[176] 13000
# nu_c[177] 12000
# nu_c[178] 9300
# nu_c[179] 9500
# nu_c[180] 20000
# nu_c[181] 20000
# nu_c[182] 10000
# nu_c[183] 20000
# nu_c[184] 6900
# nu_c[185] 20000
# nu_c[186] 12000
# nu_c[187] 12000
# nu_c[188] 9400
# nu_c[189] 7100
# nu_c[190] 7000
# nu_c[191] 12000
# nu_c[192] 9400
# nu_c[193] 12000
# nu_c[194] 20000
# nu_c[195] 12000
# nu_c[196] 6900
# nu_c[197] 14000
# nu_c[198] 12000
# nu_c[199] 7000
# nu_c[200] 12000
# nu_c[201] 9800
# nu_c[202] 7300
# nu_c[203] 20000
# nu_c[204] 7500
# nu_c[205] 9300
# nu_c[206] 16000
# nu_c[207] 14000
# nu_c[208] 12000
# nu_c[209] 9800
# nu_c[210] 10000
# nu_c[211] 11000
# nu_c[212] 12000
# nu_c[213] 6900
# nu_c[214] 11000
# nu_c[215] 7000
# nu_c[216] 20000
# nu_c[217] 11000
# nu_c[218] 7600
# nu_c[219] 20000
# nu_c[220] 20000
# nu_c[221] 7700
# nu_c[222] 20000
# nu_c[223] 9700
# nu_c[224] 13000
# nu_c[225] 7700
# nu_c[226] 7700
# nu_c[227] 12000
# nu_c[228] 9900
# nu_c[229] 20000
# nu_c[230] 8600
# nu_c[231] 6900
# nu_c[232] 20000
# nu_c[233] 9500
# nu_c[234] 10000
# nu_c[235] 6900
# nu_c[236] 20000
# nu_c[237] 9400
# nu_c[238] 9800
# nu_c[239] 8800
# nu_c[240] 6900
# nu_c[241] 11000
# nu_c[242] 12000
# nu_c[243] 9700
# nu_c[244] 9500
# nu_c[245] 11000
# nu_c[246] 15000
# nu_c[247] 9800
# nu_c[248] 9600
# nu_c[249] 9400
# nu_c[250] 7200
# nu_c[251] 9600
# nu_c[252] 12000
# nu_c[253] 12000
# nu_c[254] 12000
# nu_c[255] 12000
# nu_c[256] 11000
# nu_c[257] 12000
# nu_c[258] 11000
# nu_c[259] 11000
# nu_c[260] 7100
# nu_c[261] 11000
# nu_c[262] 7100
# nu_c[263] 13000
# nu_c[264] 20000
# nu_c[265] 17000
# nu_c[266] 9600
# nu_c[267] 12000
# nu_c[268] 7800
# nu_c[269] 12000
# nu_c[270] 20000
# nu_c[271] 9800
# nu_c[272] 7300
# nu_c[273] 11000
# nu_c[274] 20000
# nu_c[275] 12000
# nu_c[276] 20000
# nu_c[277] 12000
# nu_c[278] 9900
# nu_c[279] 6900
# nu_c[280] 14000
# nu_c[281] 12000
# nu_c[282] 12000
# nu_c[283] 11000
# nu_c[284] 9700
# nu_c[285] 11000
# nu_c[286] 9500
# nu_c[287] 7900
# nu_c[288] 6900
# nu_c[289] 9500
# nu_c[290] 9800
# nu_c[291] 16000
# nu_c[292] 9600
# nu_c[293] 7800
# nu_c[294] 9600
# nu_c[295] 9500
# nu_c[296] 14000
# nu_c[297] 7000
# nu_c[298] 14000
# nu_c[299] 20000
# nu_c[300] 7400
# nu_e[1] 20000
# nu_e[2] 20000
# nu_e[3] 20000
# nu_e[4] 20000
# nu_e[5] 20000
# nu_e[6] 20000
# nu_e[7] 20000
# nu_e[8] 20000
# nu_e[9] 20000
# nu_e[10] 20000
# nu_e[11] 20000
# nu_e[12] 20000
# nu_e[13] 20000
# nu_e[14] 20000
# nu_e[15] 20000
# nu_e[16] 20000
# nu_e[17] 20000
# nu_e[18] 20000
# nu_e[19] 20000
# nu_e[20] 20000
# nu_e[21] 20000
# nu_e[22] 20000
# nu_e[23] 20000
# nu_e[24] 20000
# nu_e[25] 20000
# nu_e[26] 20000
# nu_e[27] 20000
# nu_e[28] 20000
# nu_e[29] 20000
# nu_e[30] 20000
# nu_e[31] 20000
# nu_e[32] 20000
# nu_e[33] 20000
# nu_e[34] 20000
# nu_e[35] 20000
# nu_e[36] 20000
# nu_e[37] 20000
# nu_e[38] 20000
# nu_e[39] 20000
# nu_e[40] 20000
# nu_e[41] 20000
# nu_e[42] 20000
# nu_e[43] 20000
# nu_e[44] 20000
# nu_e[45] 20000
# nu_e[46] 20000
# nu_e[47] 20000
# nu_e[48] 20000
# nu_e[49] 20000
# nu_e[50] 20000
# nu_e[51] 20000
# nu_e[52] 20000
# nu_e[53] 20000
# nu_e[54] 20000
# nu_e[55] 20000
# nu_e[56] 20000
# nu_e[57] 20000
# nu_e[58] 20000
# nu_e[59] 20000
# nu_e[60] 20000
# nu_e[61] 20000
# nu_e[62] 20000
# nu_e[63] 20000
# nu_e[64] 20000
# nu_e[65] 20000
# nu_e[66] 20000
# nu_e[67] 20000
# nu_e[68] 20000
# nu_e[69] 20000
# nu_e[70] 20000
# nu_e[71] 20000
# nu_e[72] 20000
# nu_e[73] 20000
# nu_e[74] 20000
# nu_e[75] 20000
# nu_e[76] 20000
# nu_e[77] 20000
# nu_e[78] 20000
# nu_e[79] 20000
# nu_e[80] 20000
# nu_e[81] 20000
# nu_e[82] 20000
# nu_e[83] 20000
# nu_e[84] 20000
# nu_e[85] 20000
# nu_e[86] 20000
# nu_e[87] 20000
# nu_e[88] 20000
# nu_e[89] 20000
# nu_e[90] 20000
# nu_e[91] 20000
# nu_e[92] 20000
# nu_e[93] 20000
# nu_e[94] 20000
# nu_e[95] 20000
# nu_e[96] 20000
# nu_e[97] 20000
# nu_e[98] 20000
# nu_e[99] 20000
# nu_e[100] 20000
# nu_e[101] 20000
# nu_e[102] 20000
# nu_e[103] 20000
# nu_e[104] 20000
# nu_e[105] 20000
# nu_e[106] 20000
# nu_e[107] 20000
# nu_e[108] 20000
# nu_e[109] 20000
# nu_e[110] 20000
# nu_e[111] 20000
# nu_e[112] 20000
# nu_e[113] 20000
# nu_e[114] 20000
# nu_e[115] 20000
# nu_e[116] 20000
# nu_e[117] 20000
# nu_e[118] 20000
# nu_e[119] 20000
# nu_e[120] 20000
# nu_e[121] 20000
# nu_e[122] 20000
# nu_e[123] 20000
# nu_e[124] 20000
# nu_e[125] 20000
# nu_e[126] 20000
# nu_e[127] 20000
# nu_e[128] 20000
# nu_e[129] 20000
# nu_e[130] 20000
# nu_e[131] 20000
# nu_e[132] 20000
# nu_e[133] 20000
# nu_e[134] 20000
# nu_e[135] 20000
# nu_e[136] 20000
# nu_e[137] 20000
# nu_e[138] 20000
# nu_e[139] 20000
# nu_e[140] 20000
# nu_e[141] 20000
# nu_e[142] 20000
# nu_e[143] 20000
# nu_e[144] 20000
# nu_e[145] 20000
# nu_e[146] 20000
# nu_e[147] 20000
# nu_e[148] 20000
# nu_e[149] 20000
# nu_e[150] 20000
# nu_e[151] 20000
# nu_e[152] 20000
# nu_e[153] 20000
# nu_e[154] 20000
# nu_e[155] 20000
# nu_e[156] 20000
# nu_e[157] 20000
# nu_e[158] 20000
# nu_e[159] 20000
# nu_e[160] 20000
# nu_e[161] 20000
# nu_e[162] 20000
# nu_e[163] 20000
# nu_e[164] 20000
# nu_e[165] 20000
# nu_e[166] 20000
# nu_e[167] 20000
# nu_e[168] 20000
# nu_e[169] 20000
# nu_e[170] 20000
# nu_e[171] 20000
# nu_e[172] 20000
# nu_e[173] 20000
# nu_e[174] 20000
# nu_e[175] 20000
# nu_e[176] 20000
# nu_e[177] 20000
# nu_e[178] 20000
# nu_e[179] 20000
# nu_e[180] 20000
# nu_e[181] 20000
# nu_e[182] 20000
# nu_e[183] 20000
# nu_e[184] 20000
# nu_e[185] 20000
# nu_e[186] 20000
# nu_e[187] 20000
# nu_e[188] 20000
# nu_e[189] 20000
# nu_e[190] 20000
# nu_e[191] 20000
# nu_e[192] 20000
# nu_e[193] 20000
# nu_e[194] 20000
# nu_e[195] 20000
# nu_e[196] 20000
# nu_e[197] 20000
# nu_e[198] 20000
# nu_e[199] 20000
# nu_e[200] 20000
# nu_e[201] 20000
# nu_e[202] 20000
# nu_e[203] 20000
# nu_e[204] 20000
# nu_e[205] 20000
# nu_e[206] 20000
# nu_e[207] 20000
# nu_e[208] 20000
# nu_e[209] 20000
# nu_e[210] 20000
# nu_e[211] 20000
# nu_e[212] 20000
# nu_e[213] 20000
# nu_e[214] 20000
# nu_e[215] 20000
# nu_e[216] 20000
# nu_e[217] 20000
# nu_e[218] 20000
# nu_e[219] 20000
# nu_e[220] 20000
# nu_e[221] 20000
# nu_e[222] 20000
# nu_e[223] 20000
# nu_e[224] 20000
# nu_e[225] 20000
# nu_e[226] 20000
# nu_e[227] 20000
# nu_e[228] 20000
# nu_e[229] 20000
# nu_e[230] 20000
# nu_e[231] 20000
# nu_e[232] 20000
# nu_e[233] 20000
# nu_e[234] 20000
# nu_e[235] 20000
# nu_e[236] 20000
# nu_e[237] 20000
# nu_e[238] 20000
# nu_e[239] 20000
# nu_e[240] 20000
# nu_e[241] 20000
# nu_e[242] 20000
# nu_e[243] 20000
# nu_e[244] 20000
# nu_e[245] 20000
# nu_e[246] 20000
# nu_e[247] 20000
# nu_e[248] 20000
# nu_e[249] 20000
# nu_e[250] 20000
# nu_e[251] 20000
# nu_e[252] 20000
# nu_e[253] 20000
# nu_e[254] 20000
# nu_e[255] 20000
# nu_e[256] 20000
# nu_e[257] 20000
# nu_e[258] 20000
# nu_e[259] 20000
# nu_e[260] 20000
# nu_e[261] 20000
# nu_e[262] 20000
# nu_e[263] 20000
# nu_e[264] 20000
# nu_e[265] 20000
# nu_e[266] 20000
# nu_e[267] 20000
# nu_e[268] 20000
# nu_e[269] 20000
# nu_e[270] 20000
# nu_e[271] 20000
# nu_e[272] 20000
# nu_e[273] 20000
# nu_e[274] 20000
# nu_e[275] 20000
# nu_e[276] 20000
# nu_e[277] 20000
# nu_e[278] 20000
# nu_e[279] 20000
# nu_e[280] 20000
# nu_e[281] 20000
# nu_e[282] 20000
# nu_e[283] 20000
# nu_e[284] 20000
# nu_e[285] 20000
# nu_e[286] 20000
# nu_e[287] 20000
# nu_e[288] 20000
# nu_e[289] 20000
# nu_e[290] 20000
# nu_e[291] 20000
# nu_e[292] 20000
# nu_e[293] 20000
# nu_e[294] 20000
# nu_e[295] 20000
# nu_e[296] 20000
# nu_e[297] 20000
# nu_e[298] 20000
# nu_e[299] 20000
# nu_e[300] 20000
# s_c 20000
# s_e 15000
# deviance 20000
#
# For each parameter, n.eff is a crude measure of effective sample size,
# and Rhat is the potential scale reduction factor (at convergence, Rhat=1).
#
# DIC info (using the rule, pD = var(deviance)/2)
# pD = 7.1 and DIC = 4807.2
# DIC is an estimate of expected predictive error (lower deviance is better).The print function allows to see key posterior summaries for all parameters saved from the model, including values for posterior mean estimates, different quantiles of the posterior distribution for each parameter and diagnostic statistics such as potential scale-reduction factor or Rhat and the number of effective sample size or n.eff. Here I will not go into details about these quantities but it is enough to say that they can be used to check whether some problems occurred in the algorithm. From a first look everything seems ok. Additional checks should also be done to ensure the model behaves somewhat reasonably (e.g. no incorrect prior specification), such as posterior predictive checks, but given the simplicity of the setting I will not go into that now.
At this point however you should be asking, but what about the quantities I want to estimate, i.e. the mean QALYs and Total costs per treatment arm? how can I obtain these? Well a possible way to retrieve these is to post-process the results of the model. In particular, we can use our estimates for the conditional means \(\nu_e, \nu_c\) and standard deviations \(s_e,s_c\) in order to generate, through simulation methods, estimates for the marginal means \(\mu_e,\mu_c\) we are looking for. Although this process may seem quite complicated it is relatively simple to implement and, most importantly, can be done for most of the models we will fit, even the most complicated ones. So, in R, we do this by typing.
#obtain estimates of means by arm
#extract estimates for each mean parameter by trt group
nu_e0 <- jmodel_bn$BUGSoutput$sims.list$nu_e[,trt==0]
nu_e1 <- jmodel_bn$BUGSoutput$sims.list$nu_e[,trt==1]
nu_c0 <- jmodel_bn$BUGSoutput$sims.list$nu_c[,trt==0]
nu_c1 <- jmodel_bn$BUGSoutput$sims.list$nu_c[,trt==1]
#extract estimates for std
s_e <- jmodel_bn$BUGSoutput$sims.list$s_e
s_c <- jmodel_bn$BUGSoutput$sims.list$s_c
#create empty vectors to contain results for means by trt group
mu_e0 <- mu_c0 <- c()
mu_e1 <- mu_c1 <- c()
#set number of replications
L <- 5000
set.seed(2345) #set seed for reproducibility
#generate replications and take mean at each iteration of the posterior
for(i in 1:n.iter){
mu_e0[i] <- mean(rnorm(L,nu_e0[i,],s_e[i]))
mu_e1[i] <- mean(rnorm(L,nu_e1[i,],s_e[i]))
mu_c0[i] <- mean(rnorm(L,nu_c0[i,],s_c[i]))
mu_c1[i] <- mean(rnorm(L,nu_c1[i,],s_c[i]))
}At this point we obtained the final posterior estimates for our desired quantities, namely \(\mu_{ct}\) and \(\mu_{et}\), which can be summarised as usual, or we can even compute the incremental quantities \(\Delta_e=\mu_{e1}-\mu_{e0}\) and \(\Delta_c=\mu_{c1}-\mu_{c0}\) to see the distribution of the differences between mean outcomes by trt group (i.e. we look at the usual Cost-Effectivenss Plane).
#compute differences by arm
Delta_e <- mu_e1 - mu_e0
Delta_c <- mu_c1 - mu_c0
#plot the differences against each other
data_delta_ec <- data.frame(Delta_e,Delta_c)
ggplot(data_delta_ec, aes(x=Delta_e, y=Delta_c)) +
geom_point(size=2, shape=16) + theme_classic()
We can then produce all standard CEA output, e.g. CEAC or CE Plane, by post-processing these posterior distributions. If you want to skip the fun, we can take advantage of the R package BCEA which is dedicated to post-processing the results from a Bayesian CEA model.
#load package and provide means e and c by group as input
library(BCEA)
mu_e <- cbind(mu_e0,mu_e1)
mu_c <- cbind(mu_c0,mu_c1)
#produce CEA output
cea_res <- bcea(eff = mu_e, cost = mu_c, ref = 2)
#CE Plane (set wtp value)
ceplane.plot(cea_res, graph = "ggplot2", wtp = 10000)
#CEAC
ceac.plot(cea_res, graph = "ggplot2")# Warning: Using linewidth for a discrete variable is not advised.
#other output
summary(cea_res)#
# Cost-effectiveness analysis summary
#
# Reference intervention: intervention 2
# Comparator intervention: intervention 1
#
# Optimal decision: choose intervention 1 for k < 4900 and intervention 2 for k >= 4900
#
#
# Analysis for willingness to pay parameter k = 25000
#
# Expected net benefit
# intervention 1 15889
# intervention 2 18223
#
# EIB CEAC ICER
# intervention 2 vs intervention 1 2334.1 0.9996 4827.2
#
# Optimal intervention (max expected net benefit) for k = 25000: intervention 2
#
# EVPI 0.088877So, what you think? pretty cool…. Today we only scratch the surface of fitting Bayesian models for CEA with a very simple example based on normal distributions. In next posts I will show how these models can be tailored in a way to handle all problems of CEA data without the need to become crazy to figure out a way to fit the model or how to quantify the impact of uncertainty on the CEA results.I hope that I was able to catch your attention!