Individual level data
Want to share your content on Rbloggers? click here if you have a blog, or here if you don't.
Costeffectiveness analysis with individual level data
You can do this practical either directly in BUGS, or running the scripts from R. If you are using the BUGS graphical interface, then you will need to open the relevant model files (e.g. normalmod.txt}, cgegmod.txt), together with the relevant data (e.g. costdata.txt or costutildata.txt). Notice that for these models it is important to specify suitable initial values (particularly for the GammaGamma specification), or else BUGS may struggle to start the simulations. Conversely, if you are working directly in R, you can follow the script IPD_analysis.R}, which will guide you through the several steps and analyses required.

Analyse the cost data for the two arms of the trial. The relevant BUGS model is coded in the file normalmod.txt. Make sure you understand the code and then run it using the data stored in costdata.txt. Costs are in £1000s, so that OpenBUGS does not “overflow” during calculations with large numbers. The initial values are stored in the files normalinits1.txt, normalinits2.txt, normalinits3.txt (if you want to use 2 chains, you can select either two of these three).

Monitor the relevant nodes (decide based on your understanding of the model and the questions below. The R scripts can give you some hints if you are undecided).

What is the mean cost in each arm and the mean difference in costs between arms? (The R script also has a few more commands to manipulate the output of your simulations and produce more sophisticated and useful analyses).


In the file normalmod.txt, some code is included to calculate the deviance, a measure of model fit. The quantity 2\(\times\)log likelihood of a normal distribution is written out in BUGS model code.

Rerun the model, this time making sure you monitor the relevant nodes associated with the deviance analysis (you can skip this part if you have already monitored all the relevant nodes).

What is the deviance of the model, using this method?

Find the deviance by setting a monitor for the variable deviance.

Find the deviance using the DIC feature from the menus in OpenBUGS. .
Make sure the answers from methods b, c and d match.


Now analyse the cost and utility data for the two arms of the trial. The relevant BUGS model is coded in the file cgegmod.txt}. Inspect it and check the BUGS code makes sense.

Using the data in costutildata.txt and the initial values in the files cgeginits1.txt, cgeginits2.txt, cgeginits3.txt (if you want to use 2 chains, you can select either two of these three), run the model. If you are using BUGS, update the model for 1000 iterations. Monitor at least the nodes delta.e, delta.c, INB and CEAC. Then update the model for a further 3000 iterations.

Plot the costeffectiveness plane. .

What are the point and 95% interval estimates for the Incremental Net Benefit (INB), when the willingnessto pay is £500 per day saved in hospital? Is this consistent with the values of delta.e and delta.c _(Hint: Find the appropriate j from K[j] = (j1)*K.step, in your model code. The R file has more sophisticated code to do this)_.

What is the probability that the intensive treatment (arm 2) is costeffective for a willingnesstopay of £500 per day saved in hospital? \textit{(Hint: the R script shows all the necessary steps to compute numerically this probability. In the BUGS interface, you could modify the model code to get this number, or you could simply use the density of the relevant INB to assess this visually)._

Draw a costeffectiveness acceptability curve. (Hint: Inference \(\Rightarrow\) Compare, enter CEAC for node and K for axis, and click on scatter. Rightbuttonclick on the graphics \(\Rightarrow\) Properties. Explore the options to get a nice CEAC. The R script does this with more sophistication).

Rbloggers.com offers daily email updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/datascience job.
Want to share your content on Rbloggers? click here if you have a blog, or here if you don't.