[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: st: Stata and R

From   Stas Kolenikov <>
Subject   Re: st: Stata and R
Date   Thu, 27 Aug 2009 11:47:43 -0500

What happens to ProbCE object that you create here? That might be what
you tried to rbind, rather than INB in your last line of the cycle. I
would assume this is what you want to do.

As always, there are multiple ways of getting what you need there. You
were given a Mata solution, which looks an overkill to me, and hardly
is directly usable in further analysis. If you want a Stata variable
as a result, there are two ways of going there.

1. A somewhat clumsy way of trying to have everyting in the same data
set: assume that you have -eff- and -cost- variables for each
observation. Then you could do this:

capture noisily set obs 1000
* might give you an error message if you already have more than 1000
* we might want to capture that error message
gen INBsign = .
gen ceas = .
* this is what we will fill in

gen n = _n in 1/1000
* _n is the notation for the current observation number
* the range restriction limits the values to be from 1 to 1000

forvalues n=1/1000 {
   * n is what is known as a local macro in Stata jargon;
   * has a special reference format with opening and closing quotes
   replace INBsign = (eff*`n' - cost>0)
   sum INBsign, meanonly
   replace ceas = r(mean) in `n'
   * uses the return value from -summarize- command
   * replaces a single value in `n'-th line of data

2. What I think is a more elegant solution is to create a separate
Stata data set. In your R code, you create a separate ceas object;
unfortunately Stata can only think of one data set in memory:

postfile topost n probce using ceas, replace
* declares a new file ceas.dta to be created that is to contain two
variables, n and probce
count if ~missing(eff)
local denom = r(N)
* again, it is using the return value of the previous command
forvalues n=1/1000 {
   count if eff*`n'-cost > 0
   post topost (`n') ( r(N)/`denom' )
   * writes a new line to the ceas data set
postclose topost
use ceas, clear

HTH, Stas

On Thu, Aug 27, 2009 at 12:44 AM, Frank Peter<> wrote:
> Dear All,
> This is to thank everybody for all their answers, I have really learn a lot from the answers.
> I was trying to generate cost-effectiveness acceptability curve. The R codes I was referring to were:
> ========================================
> ceac = c()
> denom  = length(effectiveness)
> for  (n in 1000) {
>  INB = effectiveness*n - cost
>  ProbCE = length(subset(INB, INB>=0))/denom
>  ceac = rbind(ceac, c(n, INB))
> }
> ========================================
> Thanks
> Frank
> reference:
> Fenwick E, Marshall DA, Levy AR, Nichol G.     Using and interpreting cost-effectiveness acceptability curves: an example using data from a trial of management strategies for atrial fibrillation. BMC Health Serv Res. 2006 Apr 19;6:52. (

Stas Kolenikov, also found at
Small print: I use this email account for mailing lists only.

*   For searches and help try:

© Copyright 1996–2015 StataCorp LP   |   Terms of use   |   Privacy   |   Contact us   |   What's new   |   Site index