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

st: RE:comparing two densities


From   "Lachenbruch, Peter" <lachenbruch@cber.FDA.gov>
To   "statalist (E-mail)" <statalist@hsphsun2.harvard.edu>
Subject   st: RE:comparing two densities
Date   Tue, 8 Apr 2003 09:08:14 -0400

I have been facing the same problem.  There is an option that allows you to
specify the values at which you will compute the density.  
The problem I was facing was to conduct an equivalence evaluation by
comparing the density functions   I have deleted much material from the
program below in order to make the code a little simpler.  This was for a 3
group problem in which I was testing that the groups weren't equal versus an
alternative of equality.


		/* mn is the minimum value in the three groups, mx is the
maximum.  */
		scalar mn=mn1*(mn1<=mn2 & mn1<=mn3)+mn2*(mn2<=mn1 &
mn2<=mn3)+mn3*(mn3<=mn1 & mn3<=mn2)
		scalar mx=mx1*(mx1>=mx2 & mx1>=mx3)+mx2*(mx2>=mx1 &
mx2>=mx3)+mx3*(mx3>=mx1 & mx3>=mx2)
/* 	Scale so that we integrate from less than the minimum to greater
than the maximum - this allows for negative or positive minimum or maximum
*/
	scalar mn=1.25*mn*(mn<0)+.5*mn*(mn>0)
	scalar mx=1.25*mx*(mx>0)+.5*mx*(mx<0)
/* 	mn is the minimum value of the observations - choose h to be the
increment for the estimation */
	scalar h=(mx-mn)/(_N-1)
	gen xx=mn+h*(_n-1)
/* Now find densities, normalize, and find area under the maximum */
	kdensity y1,gen(xyz1 f1) at(xx) nograph
	kdensity y2,gen(xyz2 f2) at(xx) nograph
	kdensity y3,gen(xyz3 f3) at(xx) nograph
	integ f1 xx,gen(Sf1) replace
/* this is to ensure that the densities integrate to 1 - it's a small point
*/
	replace f1=f1/Sf1[_N]
	integ f2 xx,gen(Sf2) replace
	replace f2=f2/Sf2[_N]
	integ f3 xx,gen(Sf3) replace
	replace f3=f3/Sf3[_N]
/* finds maximum density */
	egen fmax=rmax(f1 f2 f3)
	integ fmax xx,gen(Sfmx) replace
/* this is compared to some critical value - I have used normal densities to
get some critical values.  However, I think it will be better to find the
critical values using a bootstrap method.  For a null of equality, it should
be much easier  */ 

Peter A. Lachenbruch, Ph. D.
Director, Division of Biostatistics
OBE/CBER/FDA
Phone    (301) 827-3320
FAX        (301) 827-5218

*
*   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