Notice: On April 23, 2014, Statalist moved from an email list to a forum, based at statalist.org.
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: st: capturing scalars in simulation
From
Nick Cox <[email protected]>
To
[email protected]
Subject
Re: st: capturing scalars in simulation
Date
Sat, 23 Feb 2013 09:51:54 +0000
In your program the first scalar you pick up from -concord- is
r(LOA_ll) and you chose to call it r(ll). Similarly in other (although
not all) cases you gave scalars different names from those used by
-concord-.
You can do that, but then you are committed to using the names you
gave elsewhere, not the names used by -concord-. This is because your
program as an r-class program overwrites previous r-class results.
it would be simpler just to run -concord- within your program and let
its r-class results be those you pick up. Otherwise what you are doing
is taking something out of one bag and putting it in another bag,
which is not needed. That way round, -mysim1- is not rclass.
-simulate- would just be picking up the results emitted by -concord-
without the need for re-bagging.
A quite different comment is that within -mysim1- it seems that a lot
of data management is done again and again, each time round the
simulation. I realise that you are taking a different sample each
time, but the rest of the management looks as if it could and should
be taken outside the program, not least because doing this for real
would, as you realise, require many more than 10 repetitions.
capture program drop mysim1
*! mysim1 Bland-Altman LOA from random pairs RJH 23021013
program define mysim1
version 12.1
drop _all
use "file location",
** data in long format with id meth replicate
*obtain random replicate for each of methods 2 & 3 (meth 1 single value)
sort id meth
by id meth : sample 1,count
** obtain requred estimates using concord program
quietly reshape wide hb repl, i(id) j(meth)
rename hb1 meth1
rename hb2 meth2
rename hb3 meth3
keep id meth*
concord meth2 meth1,
end
simulate ll =r(LOA_ll) ul =r(LOA_ul) diff =r(diff) sd =r(sd_diff) n
=r(N) , reps(10) noisily : mysim1
On Sat, Feb 23, 2013 at 8:48 AM, Richard Hiscock
<[email protected]> wrote:
> Hi all,
> I have method comparison data with haemoglobin level measured by three methods (pathology (1 reading), meth1 (3 readings) & meth2 (3 readings)) at the same time on the same subjects. Replicates are considered exchangeable within subject/method.
>
> As part of the analysis I am trying to create a dataset with a single reading for meth1 & meth2 which have been randomly sampled from within replicates for the subject & method.
> Then using concord[SJ 10-4 by T Steichen & N Cox] to obtain scalars for the mean difference & its SD + upper & lower limits of agreement
>
> I then wish to simulate this so as to obtain percentile based CI for these scalars (say 1000 reps)
> I have written the following program which successfully creates the datasets and when using noisily the analysis of each dataset occurs but I am unable to capture the scalars produced by concord except for one (the mean difference [listed as r(diff)] in concord program help file).
>
> Any guidance would be appreciated.
>
> thanks Richard Hiscock
>
> **TRY RANDOM selection from replicates
>
> capture program drop mysim1
> *! mysim1 Bland-Altman LOA from random pairs RJH 23021013
> program define mysim1, rclass
> version 12.1
> drop _all
>
> use "file location",
> ** data in long format with id meth replicate
>
> *obtain random replicate for each of methods 2 & 3 (meth 1 single value)
> sort id meth
> by id meth : sample 1,count
>
> ** obtain requred estimates using concord program
> quietly reshape wide hb repl, i(id) j(meth)
> rename hb1 meth1
> rename hb2 meth2
> rename hb3 meth3
>
> keep id meth*
> concord meth2 meth1,
> return scalar ll =r(LOA_ll)
> return scalar ul =r(LOA_ul)
> return scalar diff =r(diff)
> return scalar sd =r(sd_diff)
> return scalar n =r(N)
> end
>
> simulate ll =r(LOA_ll) ul =r(LOA_ul) diff =r(diff) sd =r(sd_diff) n =r(N) , reps(10) noisily : mysim1
*
* For searches and help try:
* http://www.stata.com/help.cgi?search
* http://www.stata.com/support/faqs/resources/statalist-faq/
* http://www.ats.ucla.edu/stat/stata/