Sergiy Radyakin <[email protected]>:
Looks like you need to take out the `rho' multiplying `sigmax_sigmay'
which means you can also lose the call to -correlate-. Or you can
keep the rho but replace `sigmax_sigmay' by the product of SDs over N.
sysuse auto, clear
svyset, srs
capture program drop st_error_of_ratio
program define st_error_of_ratio, rclass
syntax varlist(min=2 max=2)
qui svy: mean `varlist'
matrix B=e(b)
matrix V=e(V)
loc a=(B[1,1])^2*V[2,2]/(B[1,2])^4 + V[1,1]/(B[1,2])^2
loc se=sqrt(`a'-2*B[1,1]/(B[1,2])^3*V[1,2])
return scalar se=`se'
di as text "Estimated ratio (SE): " B[1,1]/B[1,2] " (" `se' ")"
end
st_error_of_ratio price length
svy: ratio price / length
On Mon, Jan 12, 2009 at 5:02 PM, Sergiy Radyakin <[email protected]> wrote:
> Dear All,
>
> I need to find a standard error of Z a ratio of two random
> variables X and Y: Z=X/Y, where the means and SEs for X and Y are
> known, as well as their corr coefficient. (In my particular case X and
> Y are numbers of people that have such and such characteristics). I am
> using delta-method according to [
> http://www.math.umt.edu/patterson/549/Delta.pdf ] (see page 2(38)). I
> then use svy:ratio to check my results and they don't match. I wonder
> if I am doing something wrong, or is it any kind of precision-related
> problem (the difference is about 7%, i.e. 1.7959 vs. 1.6795), or is my
> check simply wrong and not applicable in this case.
>
> I would appreciate if someone could look into the code and let me
> know why the results are different.
>
> Thank you,
> Sergiy Radyakin
>
> Below is a do file that one can Ctrl+C/Ctrl+V to Stata's command line:
>
> **** BEGIN OF RV_RATIO.DO ****
> sysuse auto, clear
> generate byte www=1
> svyset [pw=www]
>
> capture program drop st_error_of_ratio
> program define st_error_of_ratio, rclass
> syntax varlist(min=2 max=2)
> svy: mean `varlist'
> matrix B=e(b)
> local mux=B[1,1]
> local muy=B[1,2]
> matrix V=e(V)
> local sigma2x=V[1,1]
> local sigma2y=V[2,2]
> local sigmax_sigmay=V[1,2]
> svyset
> corr `varlist' [aw`=r(wexp)']
> local rho=r(rho)
> return scalar se = sqrt((`mux')^2*`sigma2y'/(`muy')^4 + `sigma2x'/(`muy')^2 - 2*`mux'/(`muy')^3*`rho'*`sigmax_sigmay')
> end
> st_error_of_ratio price length
> display as text "Estimated SE=" as result r(se)
> svy: ratio price / length
*
* For searches and help try:
* http://www.stata.com/help.cgi?search
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/
