Stata The Stata listserver
[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]

st: RE: RE: RE: RE: RE: Using -byable- and the -by- commands correctly?


From   "Nick Cox" <n.j.cox@durham.ac.uk>
To   <statalist@hsphsun2.harvard.edu>
Subject   st: RE: RE: RE: RE: RE: Using -byable- and the -by- commands correctly?
Date   Mon, 31 Oct 2005 17:42:20 -0000

And a slightly less rough one: 

program poissontt
	version 9
	syntax varlist(min=2 max=2), Generate(str) 
	confirm new var `generate' 
	tokenize `varlist' 
	args observed expected 
	tempvar maxprob minprob 
	
	quietly { 
		gen double `maxprob' = exp(`expected') 
		
		su `observed', meanonly 
		local max = r(max) 
		local max0 = r(max) - 1 
		
		forval I = 0/`max0' { 
			replace `maxprob' = 
`maxprob' * (`expected'^`I') / exp(lnfactorial(`I')) if `observed' <= `I' 
		}
		replace `maxprob' = 1 if `observed' == 0 
		replace `maxprob' = 1e-6 if `expected' < 1e-6
		
		gen `minprob' = 0 
		forval I = 0/`max' {
			replace `minprob' = ///
`minprob' + ((exp(-`expected')) * (`expected'^`I')) / exp(lnfactorial(`I')) ///
			if `observed' <= `I' 
		}

		gen `generate' = min(`maxprob', `minprob')   
		replace `generate' = 0 if `observed' == 0 & `expected' == 0  
		replace `generate' = max(2 * `generate', 1) 
	}		
end	

Nick 
n.j.cox@durham.ac.uk 
 
> Here is a very rough sketch: 
> 
> program poissontt
> 	version 9
> 	syntax varlist(min=2 max=2), Generate(str) 
> 	confirm new var `generate' 
> 	tempvar maxprob minprob 
> 	
> 	quietly { 
> 		gen double `maxprob' = exp(`expected') 
> 		
> 		su `observed', meanonly 
> 		local max = r(max) 
> 		local max0 = r(max) - 1 
> 		
> 		forval I = 0/`max0' { 
> 			replace `maxprob' = 
> `maxprob' * (`expected'^`I') / exp(lnfactorial(`I')) if 
> `observed' <= `I' 
> 		}
> 		replace `maxprob' = 1 if `observed' == 0 
> 		replace `maxprob' = 1e-6 if `expected' < 1e-6
> 		
> 		gen `minprob' = 0 
> 		forval I = 0/`max' {
> 			replace `minprob' = ///
> `minprob' + ((exp(-`expected')) * (`expected'^`I')) / 
> exp(lnfactorial(`I')) ///
> 			if `observed' <= `I' 
> 		}
> 
> 		gen `generate' = min(`maxprob', `minprob')   
> 		replace `generate' = 0 if `observed' == 0 & 
> `expected' == 0  
> 		replace `generate' = max(2 * `generate', 1) 
> 	}		
> end	
> 
> Nick 
> n.j.cox@durham.ac.uk 
> 
> > A broader comment is that your approach of making this 
> > byable and then looping over the data does not look 
> > the best strategy to me, as calculations look as if they can 
> > and should be recast as calculations on variables. 
> > 
> 
> *
> *   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/



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