📄 symmetry.ado
字号:
*! version 6.0.8 29sep2004
program define symmetry, rclass
version 6, missing
syntax varlist(min=2 max=2) [if] [in] [fweight] /*
*/ [, CONtrib EXact MH noTable TRend cc]
marksample touse, strok
qui {
tokenize "`varlist'"
if "`trend'"!="" {
cap confirm string var `1'
if _rc==0 {
di in red "trend not allowed with string variable `1'"
exit 198
}
cap confirm string var `2'
if _rc==0 {
di in red "trend not allowed with string variable `2'"
exit 198
}
}
if "`cc'"!="" & "`trend'"=="" {
di in red "cc option only allowed with trend option"
exit 198
}
tempvar cnt
if "`weight'" != "" {
gen double `cnt' `exp'
compress `cnt'
}
else gen byte `cnt' = 1
preserve
keep if `touse'
stack `1' `2' , into(`2') clear
keep `2'
sort `2'
qui by `2': keep if _n==1
tempfile tt
save `"`tt'"'
rename `2' `1'
Mcross `"`tt'"'
sort `1' `2'
save `"`tt'"', replace
restore, pres
keep if `touse'
keep `1' `2' `cnt'
sort `1' `2'
merge `1' `2' using `"`tt'"'
drop _merge
replace `cnt'=0 if `cnt'>=.
/* ROUTINE FOR DISPLAYING CROSS TABULATION" */
if "`table'"== "" {
noi table `1' `2',c(sum `cnt') row col center f(%5.0f)
}
/*
The rest of this code uses a strange mix of temporary
and permanent variables. The only user variables we are
using are `1' and `2'. To prevent collision, we will rename
the original variables to temporary names
*/
tempvar user1 user2
rename `1' `user1'
rename `2' `user2'
tempvar scnt
sort `user1' `user2'
by `user1' `user2': gen `scnt' = sum(`cnt')
by `user1' `user2': keep if _n == _N
/*END OF CROSS TABULATION ROUTINE*/
/* SET UP TABLE*/
tempvar gro gco
sort `user1'
by `user1': gen int `gro' = 1 if _n==1
replace `gro' = sum(`gro')
sort `user2'
by `user2': gen int `gco' = 1 if _n==1
replace `gco' = sum(`gco')
sum `gro', meanonly
local r=r(max)
rename `gro' i
rename `gco' j
rename `scnt' cv
sort i j
*keep i j cv
gen int ncv = cv[_n+(`r'-1)*(j-i)]
if "`trend'"!="" {
if "`cc'"!="" {
MTrend `user1' `user2' ccor
}
else MTrend `user1' `user2'
local ptrend=chiprob(1, r(trend))
local trend=r(trend)
}
keep i j cv ncv
/* END SET UP TABLE */
/* ROUTINE FOR MARGINAL HOMOGENEITY*/
sort i
gen diag=cv if i==j
by i: gen mr=sum(cv)
by i: replace mr=mr[_N]
gen mmr=mr-diag
sort j
by j: gen mc=sum(cv)
by j: replace mc=mc[_N]
gen mmc=mc-diag
gen nd=mc-mr
sort i j
tempname d
mkmat nd if i==j,mat(`d')
drop nd
drop diag
if "`mh'"=="mh" {
gen double m = ((mr-mc)^2)/(mr+mc) if i==j
gen double mm = ((mmr-mmc)^2)/(mmr+mmc) if i==j
sum m, meanonly
local k = r(N)
local m = r(sum)
local df2 = `k'-1
gen tdf = 1 if i==j
replace tdf =. if mmr==0 & mmc==0
count if tdf==1
drop tdf
local df3 = r(N)-1
sum mm, meanonly
local mm = r(sum)
local p2 = chiprob(`df2',`m')
local p3 = chiprob(`df3',(`k'-1)*`mm'/`k')
}
/* ROUTINE FOR STUART-MAXWELL MARGINAL HOMOGENEITY*/
tempname V Vi U X
gen vii = mmr+mmc
drop mmr mmc
replace vii = -1*(cv+ncv) if i~=j
matrix `V' = J(`r',`r',0)
sort i j
local i 1
local ob 1
while `i'<=`r' {
local j 1
while `j'<=`r' {
mat `V'[`i',`j']= vii[`ob']
local ob=`ob'+1
local j=`j'+1
}
local i=`i'+1
}
mat `Vi'=syminv(`V')
mat `U'=`d''*`Vi'
mat `X'=`U'*`d'
local bmd=`X'[1,1]
local df1=`r'-1
local p4=chiprob(`df1',`bmd')
/* ROUTINE FOR SYMMETRY TEST */
gen double sym = ((cv-ncv)^2)/(-1*vii) if i<j
drop vii
replace sym=0 if sym>=. & i<j
sort i j
if "`contrib'"=="contrib" {
local i 1
local ob 1
noi di
noi di in green " Contribution"
noi di in green " to symmetry"
noi di in green " Cells chi-squared"
noi di in smcl in green "{hline 14} {hline 14}"
while `ob'<=_N {
if sym[`ob']<. {
noi di " n" i[`ob'] "_" j[`ob'] " & n" /*
*/ j[`ob'] "_" i[`ob'] " " %9.4f sym[`ob']
}
local ob=`ob'+1
}
}
gen d=0 if i<j
replace d=1 if cv==0 & ncv==0 & d==0
gen long nd=sum(d)
sum sym, meanonly
local x2 = r(sum)
local df = (`r'*(`r'-1)/2)-nd[_N]
local p = chiprob(`df',`x2')
drop d nd
/* END SYMMETRY TEST */
sum cv, meanonly
local totpair = r(sum)
/* EXACT TEST ALGORITHM nxn tables */
gen pval=.
if "`exact'"!="" {
keep if i<j
keep i j cv ncv
rename cv v1
rename ncv v2
*order i j v1 v2
sort i j
gen den= v1+v2
/* OBSERVED TABLE'S probability */
tempname obsval
gen double oval=comb(den,v1)*(0.5)^den
gen double obval=oval
replace obval=oval*obval[_n-1] if _n>1
scalar `obsval'=obval[_N]
drop obval oval
gen double pval = . /* if obsval>=. */
if `obsval'<. {
gen long nv1=den
gen long nv2=0
gen u = (den - (mod(den,2)~=0)) / 2
gen negu= -u
sort negu i j
drop negu
local i 1
while `i'<=_N {
local den`i'=den[`i']
local u`i'=u[`i'] /* upper bound per variable */
local i=`i'+1
}
local r=_N /* variables go from 1 to r=# of collumns */
drop _all /* First Observation */
set obs 1 /* set to zero vector */
local i 1
while `i'<=`r'{
gen v`i'= 0
local i = `i'+1
}
expand `u1' + 1
replace v1 = _n - 1
local i 2
local myval "v1"
gen double nnp = 1
while `i' <= `r' {
expand `u`i'' + 1
local myval "`myval' v`i'"
sort `myval'
by `myval': replace v`i'=_n-1
replace nnp=1
local j 1
while `j'<=`r' {
replace nnp= nnp /*
*/ * comb(`den`j'',v`j')*(0.5)^`den`j''
local j=`j'+1
}
drop if nnp>`obsval'
local i=`i'+1
}
drop nnp
local i 1
gen double np=1
gen eq=0
while `i'<=`r' {
replace np=np*comb(`den`i'',v`i')*(0.5)^`den`i''
replace eq=eq+1 if `den`i'' - v`i' == v`i'
local i=`i'+1
}
keep if np< = `obsval'
replace np=np*2^(`r'-eq)
gen double pval=sum(np)
}
} /*end of exact for nxn */
noi di
noi di in green _col(8) /*
*/ " chi2 df Prob>chi2"
noi di in smcl in gr "{hline 38}{c TT}{hline 33}"
noi di in smcl in gr "Symmetry (asymptotic) {c |} " /*
*/ in ye %9.2f `x2' _col(56) `df' _col(66) %5.4f `p'
noi di in smcl in gr "Marginal homogeneity (Stuart-Maxwell) {c |} " /*
*/ in ye %9.2f `bmd' /*
*/ _col(56) `df1' _col(66) %5.4f `p4'
if "`mh'"!="" {
noi di in smcl in gr /*
*/ "Marginal homogeneity (Bickenboller) {c |} " /*
*/ in ye %9.2f `m' /*
*/ _col(56) `df2' _col(66) %5.4f `p2'
noi di in smcl in gr /*
*/ "Marginal homogeneity (no diagonals) {c |} " /*
*/ in ye %9.2f (`k'-1)*`mm'/`k' /*
*/ _col(56) `df3' _col(66) %5.4f `p3'
ret scalar p_nd = `p3'
ret scalar df_nd = `df3'
ret scalar chi2_nd = (`k'-1)*`mm'/`k'
ret scalar p_b = `p2'
ret scalar df_b = `df2'
ret scalar chi2_b = `m'
}
if "`trend'" == "" { noi di in smcl in gr "{hline 38}{c BT}{hline 33}" }
else {
noi di in smcl in gr "{hline 38}{c +}{hline 33}"
noi di in smcl in gr /*
*/ "Linear trend in the (log) RR {c |} " /*
*/ in ye %9.2f `trend' /*
*/ _col(56) "1" _col(66) %5.4f `ptrend'
noi di in smcl in gr "{hline 38}{c BT}{hline 33}"
ret scalar chi2_t= `trend'
ret scalar p_trend= `ptrend'
}
if "`exact'"!="" {
noi di in gr "Symmetry (exact significance probability) " /*
*/ in ye _col(66) %5.4f pval[_N]
if pval[_N]>=. {
noi di in gr "(Too many observations for exact test," /*
*/ " asymptotic result should hold)"
}
ret scalar p_exact= pval[_N]
}
ret scalar p_sm = `p4'
ret scalar df_sm = `df1'
ret scalar chi2_sm = `bmd'
ret scalar p = `p'
ret scalar df = `df'
ret scalar chi2 = `x2'
ret scalar N_pair = `totpair'
}
end
program define Mcross /* using `tt' */
args using
local nob = _N
tempfile cross2
tempvar order midx
preserve
quietly use `"`using'"', clear
quietly {
gen long `order'=_n
expand `nob', clear
sort `order'
by `order': gen long `midx' = _n
sort `midx' `order'
drop `order'
save `"`cross2'"', replace
restore, preserve
gen long `midx' = _n
sort `midx'
merge `midx' using "`cross2'"
drop `midx' _merge
restore, not
}
end
program define MTrend, rclass /* r(trend)*/
args case control ccor
if "`ccor'"!="" {
local cc=0.5
}
else local cc 0
tempvar diff dexp N num snum den sden
qui gen long `diff' = cv-ncv if i<j
qui gen `dexp' = `case' - `control' if i<j
qui gen `N' = cv+ncv if i<j
qui gen double `num' =`diff'*`dexp'
qui gen double `snum'=sum(`num')
qui replace `num'=(`snum'[_N] - `cc')^2
qui gen double `den'=`N'*(`dexp'^2)
qui gen double `sden'=sum(`den')
ret scalar trend =`num'/`sden'[_N]
end
exit
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -