Certainly Nick,
please see below as requested. This program is based on Fortran
functions incorporated in a Person-Years analysis program.
capture program drop poissontt
program poissontt, sortpreserve byable(recall, noheader)
/* make program usable with 'by' command - see manual [P] page 12 */
version 9
quietly
local poissval = 0
local minprob = 0
local maxprob = 0
local observed = `1' /* O */
local expected = `2' /* E */
local replacevar = "replace `3' = "
/* probability Observed no >= O */
if `observed' == 0 {
local maxprob = 1
}
else if `expected' < 0.000001 {
local maxprob = 0.000001
}
else {
local I = 0
local sumcalc = 0
local tempcalc = 0
local endloop = `observed' - 1
while `I' <= `endloop' {
local tempcalc = ((exp(-1 * `expected')) *
(`expected' ^ `I')) / exp(lnfactorial(`I'))
local sumcalc = `sumcalc' + `tempcalc'
local I = `I' + 1
}
local maxprob = 1 -`sumcalc'
}
/* probability Observed no <= O */
local I = 0
local sumcalc = 0
local tempcalc = 0
local endloop = `observed'
while `I' <= `endloop' {
local tempcalc = ((exp(-1 * `expected')) * (`expected' ^
`I')) / exp(lnfactorial(`I'))
local sumcalc = `sumcalc' + `tempcalc'
local I = `I' + 1
}
local minprob = `sumcalc'
/* evaluate probsbility to calculate statistic */
if `observed' == 0 & `expected' == 0 {
local poissval = 0 /* empty */
}
else {
if `minprob' <= `maxprob' {
local poissval = `minprob'
}
else {
local poissval = `maxprob'
}
local poissval = 2 * `poissval'
if `poissval' > 1 {
local poissval = 1
}
}
`replacevar' `poissval'
noisily
end
-----Original Message-----
From: [email protected]
[mailto:[email protected]] On Behalf Of Nick Cox
Sent: 31 October 2005 10:39
To: [email protected]
Subject: st: RE: Using -byable- and the -by- commands correctly?
Please show us the whole program.
Nick
[email protected]
David Winter
> Firstly, thanks to Svend Juul for the helpful comments on survival
> analysis.
>
> The following 4 records are part of a 27 record dataset sorted by
> ageband (5,10,15,20..85) within sex (MALE,FEMALE).
>
> Sex Ageband Observed Expected PoissStat
> MALE 5 5 .2585
> MALE 10 11 .6054
> MALE 15 6 .6516
> MALE 20 5 .6435
>
> I have written a program to calculate a two-tailed poisson statistic,
> that has three arguements:
>
> `1' = Observed number
> `2' = Expected number
> `3' = variable to hold the poisson statistic
>
> I have made the program -byable- :
>
> Program poissontt, sortpreserve byable(recall, no header)
> ..
> ..
> local poissval = 0 /* macro to store the result */
> ..
> local observed = `1'
> local expected =`2'
> local replacevar = "replace `3' = "
> ..
> <define poissval>
> ..
> `replacevar' `poissval'
> end
>
> In a -.do- file, I run the commands
>
> gen PoissStat = .
> by sex ageband: poissontt Observed Expected PoissStat
>
> Stata diligently loops through the data and using -set trace- I have
> tracked the iterations. It does not appear to pick up the values of
> Observed and Expected for each record, rather those of the
> first record
> only. The calculation is correct for record 1 (Obs=5, Expect=.258,
> poisson(two-tailed=.0000153) and the PoissStat variable is updated in
> each record, but with the value for record 1 only.
>
> Is there something I am missing either in the program or in using the
> -by- command?
*
* For searches and help try:
* http://www.stata.com/support/faqs/res/findit.html
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/
*
* For searches and help try:
* http://www.stata.com/support/faqs/res/findit.html
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/
