You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
919 lines
27 KiB
Plaintext
919 lines
27 KiB
Plaintext
7 months ago
|
*! metabias version 1.2.2 07feb01 TJS (STB-61: sbe19.4)
|
||
|
|
program define metabias, rclass
|
||
|
|
version 6.0
|
||
|
|
* version 1.0.0 24oct97 TJS STB-41 sbe19
|
||
|
|
* version 1.1.2 21apr98 TJS Unweighted Egger analysis only STB-44 sbe19.1
|
||
|
|
* version 1.1.3 1oct99 TJS fixed stratified bug STB-44 sbe19.1
|
||
|
|
* version 1.2.0 17feb00 TJS Version 6 and gweight option (STB-57: sbe19.2)
|
||
|
|
* version 1.2.1 31jan01 TJS touched for version 7.0
|
||
|
|
*
|
||
|
|
*! syntax: metabias varlist [if] [in]
|
||
|
|
*! [ , Graph(str) GWeight Vari CI BY(varname) LEvel(real 95) * ]
|
||
|
|
|
||
|
|
if ("`*'" == "") {
|
||
|
|
di "Syntax is:"
|
||
|
|
di in wh "metabias " in gr "{ theta { se | var } | " _c
|
||
|
|
di in gr "exp(theta) | ll ul [cl] } [" _c
|
||
|
|
di in wh "if " in gr "exp] [" in wh "in " in gr "range]"
|
||
|
|
di in gr " [ " in wh ", by(" in gr "by_var"in wh ")" _c
|
||
|
|
di in gr " { " in wh "v" in gr "ar | " in wh "ci" in gr " } " _c
|
||
|
|
di in wh "g" in gr "raph" in wh "(" in gr "{ " in wh "b" _c
|
||
|
|
di in gr "egg | " in wh "e" in gr "gger }" in wh ")"
|
||
|
|
di in wh " gw" in gr "eight" in wh " lev" _c
|
||
|
|
di in gr "el" in wh "(" in gr "#" in wh ") " in gr "graph_options ]"
|
||
|
|
di _n in gr " where { a | b |...} means choose" _c
|
||
|
|
di in gr " one and only one of {a, b, ...}"
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
|
||
|
|
* Setup
|
||
|
|
syntax varlist(numeric min=2 max=4) [if] [in] /*
|
||
|
|
*/ [ , Graph(str) Vari CI BY(varname) LEvel(real 95) * ]
|
||
|
|
* Note: option GWeight is passed via * directly to the graphics subs
|
||
|
|
marksample touse
|
||
|
|
tokenize `varlist'
|
||
|
|
|
||
|
|
tempvar byg theta setheta var w sw vt wtheta swtheta
|
||
|
|
tempvar zz Ts wl swl RRm bylabel
|
||
|
|
tempname k ks sdks p zu pcc zcc c sks svks sk oe sbv sv
|
||
|
|
* v7 fix
|
||
|
|
tempvar bylabl
|
||
|
|
|
||
|
|
local theta `1'
|
||
|
|
if "`3'" == "" { local setheta `2' }
|
||
|
|
else {
|
||
|
|
tempvar ll ul cl
|
||
|
|
local ll `2'
|
||
|
|
local ul `3'
|
||
|
|
local cl `4'
|
||
|
|
}
|
||
|
|
|
||
|
|
* input error traps
|
||
|
|
if "`ci'" != "" & "`vari'" != "" {
|
||
|
|
di _n in re "Error: options 'ci' and 'var' cannot " _c
|
||
|
|
di in re "be specified together."
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if "`ci'" == "ci" & "`ul'" != "" {
|
||
|
|
di _n in bl "Note: option 'ci' specified."
|
||
|
|
}
|
||
|
|
if "`ci'" == "ci" & "`ul'" == "" {
|
||
|
|
di _n in re "Error: option 'ci' specified but varlist " _c
|
||
|
|
di in re "has only 2 variables."
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if "`ci'" != "ci" & "`vari'" != "vari" & "`ul'" != "" {
|
||
|
|
di _n in bl "Warning: varlist has 3 variables but option " _c
|
||
|
|
di in bl "'ci' not specified; 'ci' assumed."
|
||
|
|
local ci "ci"
|
||
|
|
local vari ""
|
||
|
|
}
|
||
|
|
if "`vari'" == "vari" & "`ul'" != "" {
|
||
|
|
di _n in re "Error: option 'var' specified but varlist " _c
|
||
|
|
di in re "has more than 2 variables."
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if "`vari'" == "vari" & "`ul'" == "" {
|
||
|
|
di _n in bl "Note: option 'var' specified."
|
||
|
|
}
|
||
|
|
if "`vari'" != "vari" & "`ul'" == "" {
|
||
|
|
di _n in bl "Note: default data input format (theta, " _c
|
||
|
|
di in bl "se_theta) assumed."
|
||
|
|
}
|
||
|
|
|
||
|
|
* Select data to analyze
|
||
|
|
if "`ul'" == "" { markout `touse' `theta' `setheta' }
|
||
|
|
else { markout `touse' `theta' `ll' `ul' }
|
||
|
|
|
||
|
|
preserve /* Note: data preserved here */
|
||
|
|
|
||
|
|
* Generate `by' groups
|
||
|
|
if "`by'" != "" {
|
||
|
|
confirm var `by'
|
||
|
|
sort `by'
|
||
|
|
qui by `by': gen byte `byg' = _n==1
|
||
|
|
qui replace `byg' = sum(`byg')
|
||
|
|
local byn = `byg'[_N]
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
qui gen byte `byg' = 1
|
||
|
|
local byn = 1
|
||
|
|
}
|
||
|
|
|
||
|
|
* Generate `by' labels -- if required
|
||
|
|
if ("`by'" != "") {
|
||
|
|
capture decode `by', gen(`bylabel')
|
||
|
|
if _rc != 0 {
|
||
|
|
local type : type `by'
|
||
|
|
if substr("`type'",1,3) != "str" {
|
||
|
|
qui gen str8 `bylabel' = string(`by')
|
||
|
|
}
|
||
|
|
else { qui gen `type' `bylabel' = `by' }
|
||
|
|
}
|
||
|
|
* v7 fix
|
||
|
|
encode `bylabel', gen(`bylabl')
|
||
|
|
}
|
||
|
|
|
||
|
|
* Do calculations
|
||
|
|
* initial calculations...
|
||
|
|
if "`vari'" == "vari" { qui replace `setheta' = sqrt(`setheta')}
|
||
|
|
|
||
|
|
if "`ci'" == "ci" {
|
||
|
|
capture confirm variable `cl'
|
||
|
|
if _rc~=0 { qui gen `zz' = invnorm(.975) }
|
||
|
|
else {
|
||
|
|
qui replace `cl' = `cl' * 100 if `cl' < 1
|
||
|
|
qui gen `zz' = -1 * invnorm((1- `cl' / 100) / 2 )
|
||
|
|
qui replace `zz' = invnorm(.025) if `zz'==.
|
||
|
|
}
|
||
|
|
qui gen `setheta' = ( ln(`ul') - ln(`ll')) / 2 / `zz'
|
||
|
|
qui replace `theta' = ln(`theta')
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`if'" != "" { ifexp "`if'" }
|
||
|
|
|
||
|
|
if "`by'" != "" {
|
||
|
|
scalar `sk' = 0
|
||
|
|
scalar `sks' = 0
|
||
|
|
scalar `svks' = 0
|
||
|
|
scalar `sbv' = 0
|
||
|
|
scalar `sv' = 0
|
||
|
|
}
|
||
|
|
|
||
|
|
* loop through by-values
|
||
|
|
local j = 1
|
||
|
|
while `j' <= `byn' { /* start of loop for each `by' group */
|
||
|
|
|
||
|
|
summ `touse' if `touse' & `byg' == `j', meanonly
|
||
|
|
local data = _result(1)
|
||
|
|
|
||
|
|
* Calculate stats
|
||
|
|
qui {
|
||
|
|
gen `var' = `setheta'^2
|
||
|
|
gen `w' = 1/`var'
|
||
|
|
egen `sw' = sum(`w') if `touse' & `byg' == `j'
|
||
|
|
gen `vt' = `var' - 1 / `sw'
|
||
|
|
gen `wtheta' = `w' * `theta'
|
||
|
|
egen `swtheta' = sum(`wtheta') if `touse' & `byg' == `j'
|
||
|
|
gen `Ts' = (`theta' - `swtheta' / `sw') / sqrt(`vt')
|
||
|
|
gen `wl' = `w' * `theta'
|
||
|
|
egen `swl' = sum(`wl') if `touse' & `byg' == `j'
|
||
|
|
gen `RRm' = `swl' / `sw'
|
||
|
|
scalar `oe' = `RRm'
|
||
|
|
}
|
||
|
|
|
||
|
|
qui capture ktau2 `var' `Ts' if `touse' & `byg' == `j'
|
||
|
|
if _rc == 0 {
|
||
|
|
scalar `k' = $S_1
|
||
|
|
scalar `ks' = $S_4
|
||
|
|
scalar `sdks' = $S_5
|
||
|
|
scalar `p' = $S_6
|
||
|
|
scalar `zu' = $S_7
|
||
|
|
scalar `pcc' = $S_8
|
||
|
|
scalar `zcc' = $S_9
|
||
|
|
scalar `c' = $S_10
|
||
|
|
}
|
||
|
|
else if _rc == 2001 {
|
||
|
|
scalar `k' = `data'
|
||
|
|
scalar `ks' = .
|
||
|
|
scalar `sdks' = .
|
||
|
|
scalar `p' = .
|
||
|
|
scalar `zu' = .
|
||
|
|
scalar `pcc' = .
|
||
|
|
scalar `zcc' = .
|
||
|
|
scalar `c' = .
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
di in re "error " _rc " in call to ktau2"
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`by'" != "" & `k' > 1 {
|
||
|
|
scalar `sk' = `sk' + `k'
|
||
|
|
scalar `sks' = `sks' + `ks'
|
||
|
|
scalar `svks' = `svks' + `sdks'^2
|
||
|
|
}
|
||
|
|
|
||
|
|
* Egger's bias test
|
||
|
|
tempvar prec snd
|
||
|
|
qui gen `prec'= 1 / `setheta'
|
||
|
|
qui gen `snd' = `theta' / `setheta'
|
||
|
|
qui regr `snd' `prec' if `touse' & `byg' == `j'
|
||
|
|
capture matrix b = get(_b)
|
||
|
|
if _rc == 0 {
|
||
|
|
local df = e(N) - 2
|
||
|
|
local bias = b[1,2]
|
||
|
|
capture matrix V = get(VCE)
|
||
|
|
if _rc == 0 {
|
||
|
|
local pb = tprob(`df', b[1,2] / sqrt(V[2,2]))
|
||
|
|
* 1.1.3 fix
|
||
|
|
* if "`by'" != "" {
|
||
|
|
if "`by'" != "" & V[2,2] != 0 & `data' > 0 {
|
||
|
|
scalar `sbv' = `sbv' + `bias' / V[2,2]
|
||
|
|
scalar `sv' = `sv' + 1 / V[2,2]
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
local bias = .
|
||
|
|
local pb = .
|
||
|
|
}
|
||
|
|
|
||
|
|
* Display results
|
||
|
|
if "`by'" != "" {
|
||
|
|
* use this display if a by_var was specified...
|
||
|
|
* .....display output header
|
||
|
|
if `j' == 1 {
|
||
|
|
di " "
|
||
|
|
di in gr "Tests for Publication Bias"
|
||
|
|
di " "
|
||
|
|
local sp = 8 - length("`by'")
|
||
|
|
#delimit ;
|
||
|
|
di in gr "-------------------------------------------------"
|
||
|
|
"------------------------------" ;
|
||
|
|
di in gr " | | Begg's Begg's"
|
||
|
|
" cont. corr. | Egger's " ;
|
||
|
|
di in gr _skip(`sp') "`by' | n | score s.d. z"
|
||
|
|
" p z p | bias p" ;
|
||
|
|
di in gr "---------+------+--------------------------------"
|
||
|
|
"---------------+--------------" ;
|
||
|
|
#delimit cr
|
||
|
|
local scs " "
|
||
|
|
}
|
||
|
|
* v7 fix
|
||
|
|
* local sp = 8 - length(`bylabel')
|
||
|
|
local blab: label (`bylabl') `j'
|
||
|
|
local sp = 8 - length("`blab'")
|
||
|
|
|
||
|
|
if `c' == 1 {
|
||
|
|
local cs "*"
|
||
|
|
local scs "*"
|
||
|
|
}
|
||
|
|
else {local cs " "}
|
||
|
|
* .....display results for each by-value
|
||
|
|
if `data' > 0 {
|
||
|
|
* v7 fix
|
||
|
|
di in gr _skip(`sp') "`blab' | " in ye %4.0f `k' _c
|
||
|
|
* di in gr _skip(`sp') `bylabel' " | " in ye %4.0f `k' _c
|
||
|
|
di in gr " |" in ye %4.0f `ks' in gr "`cs' " _c
|
||
|
|
di in ye %6.3f `sdks' " " %6.2f `zu' " " %6.3f `p' " " _c
|
||
|
|
di in ye %6.2f `zcc' " " %6.3f `pcc' in gr " |" _c
|
||
|
|
di in ye %6.2f `bias' " " %6.3f `pb'
|
||
|
|
}
|
||
|
|
* .....do stratified calculations
|
||
|
|
if `j' == `byn' {
|
||
|
|
scalar `zu' = `sks' / sqrt(`svks')
|
||
|
|
scalar `p' = 2 * (1 - normprob(abs(`zu')))
|
||
|
|
scalar `zcc' = sign(`sks')*(abs(`sks') - 1) / sqrt(`svks')
|
||
|
|
scalar `pcc' = 2 * (1 - normprob(abs(`zcc')))
|
||
|
|
scalar `k' = `sk'
|
||
|
|
scalar `ks' = `sks'
|
||
|
|
scalar `sdks' = sqrt(`svks')
|
||
|
|
drop `sw' `wl' `swl' `RRm'
|
||
|
|
qui egen `sw' = sum(`w') if `touse'
|
||
|
|
qui gen `wl' = `w' * `theta'
|
||
|
|
qui egen `swl' = sum(`wl') if `touse'
|
||
|
|
qui gen `RRm' = `swl' / `sw'
|
||
|
|
scalar `oe' = `RRm'
|
||
|
|
|
||
|
|
qui regr `snd' `prec' if `touse'
|
||
|
|
capture matrix b = get(_b)
|
||
|
|
if _rc == 0 {
|
||
|
|
matrix V = get(VCE)
|
||
|
|
local df = e(N) - 2
|
||
|
|
local bias = b[1,2]
|
||
|
|
local pb = tprob(`df', b[1,2]/sqrt(V[2,2]))
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
local bias = .
|
||
|
|
local pb = .
|
||
|
|
}
|
||
|
|
|
||
|
|
local bias = `sbv' / `sv'
|
||
|
|
local pb = 1 - normprob(`bias' / sqrt(1 / `sv'))
|
||
|
|
|
||
|
|
* .....and display overall (stratified) results
|
||
|
|
di in gr "---------+------+----------------------------" _c
|
||
|
|
di in gr "-------------------+--------------"
|
||
|
|
di in gr " overall | " in ye %4.0f `sk' _c
|
||
|
|
di in gr " |" in ye %4.0f `sks' in gr "`scs' " _c
|
||
|
|
di in ye %6.3f sqrt(`svks') " " %6.2f `zu' " " %6.3f `p' _c
|
||
|
|
di in ye " " %6.2f `zcc' " " %6.3f `pcc' in gr " |" _c
|
||
|
|
di in ye %6.2f `bias' " " %6.3f `pb'
|
||
|
|
di in gr "---------------------------------------------" _c
|
||
|
|
di in gr "----------------------------------"
|
||
|
|
if "`scs'" == "*" {
|
||
|
|
di in gr _skip(21) "`scs' (corrected for ties)"
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
else {
|
||
|
|
|
||
|
|
* use this display if no by_var was specified...
|
||
|
|
* Begg's
|
||
|
|
di _n in gr "Tests for Publication Bias"
|
||
|
|
di _n in gr "Begg's Test"
|
||
|
|
di " "
|
||
|
|
di in gr " adj. Kendall's Score (P-Q) = " in ye %7.0f `ks'
|
||
|
|
di _c in gr " Std. Dev. of Score = " in ye %7.2f `sdks'
|
||
|
|
if `c' == 1 { di in gr " (corrected for ties)" }
|
||
|
|
else { di " " }
|
||
|
|
di in gr " Number of Studies = " in ye %7.0f `k'
|
||
|
|
di in gr " z = " in ye %7.2f `zu'
|
||
|
|
di in gr " Pr > |z| = " in ye %7.3f `p'
|
||
|
|
di _c in gr " z = " in ye %7.2f `zcc'
|
||
|
|
di in gr " (continuity corrected)"
|
||
|
|
di _c in gr " Pr > |z| = " in ye %7.3f `pcc'
|
||
|
|
di in gr " (continuity corrected)"
|
||
|
|
|
||
|
|
* Egger's
|
||
|
|
tempvar prec snd
|
||
|
|
qui gen `prec'= 1 / `setheta'
|
||
|
|
qui gen `snd' = `theta' / `setheta'
|
||
|
|
qui regr `snd' `prec' if `touse'
|
||
|
|
capture matrix b = get(_b)
|
||
|
|
if _rc == 0 {
|
||
|
|
matrix V = get(VCE)
|
||
|
|
local obs = e(N)
|
||
|
|
local df = `obs' - 2
|
||
|
|
matrix colnames b = slope bias
|
||
|
|
matrix rownames V = slope bias
|
||
|
|
matrix colnames V = slope bias
|
||
|
|
matrix post b V, dep(Std_Eff) dof(`df') obs(`obs')
|
||
|
|
di _n in gr "Egger's test"
|
||
|
|
matrix mlout, level(`level')
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
di _n in gr "Egger's Test" _n
|
||
|
|
di in ye " - undefined for only 1 study"
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
cap drop `var' `w' `sw' `vt' `wtheta' `swtheta'
|
||
|
|
cap drop `Ts' `wl' `swl' `RRm'
|
||
|
|
local j = `j' + 1
|
||
|
|
|
||
|
|
} /* end of loop for each `by' group */
|
||
|
|
|
||
|
|
* Graph a bias plot
|
||
|
|
local g = lower(substr("`graph'",1,1))
|
||
|
|
if "`g'" == "b" {
|
||
|
|
beggph `theta' `setheta' `touse', level(`level') `ci' `options'
|
||
|
|
}
|
||
|
|
if "`g'" == "e" {
|
||
|
|
egggph `theta' `setheta' `touse', level(`level') `options'
|
||
|
|
}
|
||
|
|
|
||
|
|
* Save globals
|
||
|
|
global S_1 = `k'
|
||
|
|
global S_2 = `ks'
|
||
|
|
global S_3 = `sdks'
|
||
|
|
global S_4 = `p'
|
||
|
|
global S_5 = `pcc'
|
||
|
|
global S_6 = `bias'
|
||
|
|
global S_7 = `pb'
|
||
|
|
global S_8 = `oe'
|
||
|
|
|
||
|
|
* return globals
|
||
|
|
return scalar k = `k'
|
||
|
|
return scalar score = `ks'
|
||
|
|
return scalar score_sd = `sdks'
|
||
|
|
return scalar Begg_p = `p'
|
||
|
|
return scalar Begg_pcc = `pcc'
|
||
|
|
return scalar Egger_bc = `bias'
|
||
|
|
return scalar Egger_p = `pb'
|
||
|
|
return scalar effect = `oe'
|
||
|
|
|
||
|
|
exit
|
||
|
|
end
|
||
|
|
|
||
|
|
* ***************************************************
|
||
|
|
|
||
|
|
program define beggph
|
||
|
|
version 6.0
|
||
|
|
|
||
|
|
* creates the Begg funnel plot graph
|
||
|
|
|
||
|
|
* Setup
|
||
|
|
|
||
|
|
syntax varlist(min=3 max=3) [if] [in] [, CI L1title(str) /*
|
||
|
|
*/ L2title(str) Connect(str) Symbol(str) SOrt Pen(str) /*
|
||
|
|
*/ T2title(str) B2title(str) YLAbel(str) XLAbel(str) /*
|
||
|
|
*/ LEVel(integer $S_level) GAp(str) GWeight * ]
|
||
|
|
tokenize `varlist'
|
||
|
|
|
||
|
|
tempvar touse theta setheta
|
||
|
|
|
||
|
|
local theta `1'
|
||
|
|
local setheta `2'
|
||
|
|
local touse `3'
|
||
|
|
|
||
|
|
preserve
|
||
|
|
|
||
|
|
* Graph options
|
||
|
|
if "`connect'" == "" { local connect "lll." }
|
||
|
|
else {
|
||
|
|
local lll = length("`connect'")
|
||
|
|
if `lll' == 1 { local connect "`connect'll." }
|
||
|
|
else if `lll' == 2 { local connect "`connect'l." }
|
||
|
|
else if `lll' == 3 { local connect "`connect'." }
|
||
|
|
}
|
||
|
|
local connect "co(`connect')"
|
||
|
|
|
||
|
|
if "`symbol'" == "" { local symbol "iiio" }
|
||
|
|
else {
|
||
|
|
local lll = length("`symbol'")
|
||
|
|
if `lll' == 1 { local symbol "`symbol'iio" }
|
||
|
|
else if `lll' == 2 { local symbol "`symbol'io" }
|
||
|
|
else if `lll' == 3 { local symbol "`symbol'o" }
|
||
|
|
}
|
||
|
|
local symbol "sy(`symbol')"
|
||
|
|
|
||
|
|
if "`pen'" == "" { local pen "3552" }
|
||
|
|
else {
|
||
|
|
local lll = length("`pen'")
|
||
|
|
if `lll' == 1 { local pen "`pen'552" }
|
||
|
|
else if `lll' == 2 {
|
||
|
|
local pen = "`pen'" + substr("`pen'",2,1) + "2"
|
||
|
|
}
|
||
|
|
else if `lll' == 3 { local pen "`pen'2" }
|
||
|
|
}
|
||
|
|
local pen "pen(`pen')"
|
||
|
|
|
||
|
|
if "`l2title'" == "" {
|
||
|
|
local l2title : variable label `theta'
|
||
|
|
if "`l2title'" == "" { local l2title "`theta'" }
|
||
|
|
}
|
||
|
|
if "`ci'" == "" { local l2title "l2(`l2title')" }
|
||
|
|
else { local l2title "l2(log[`l2title'])" }
|
||
|
|
|
||
|
|
if "`l1title'" == "" { local l1title "" "" }
|
||
|
|
local l1title "l1(`l1title')"
|
||
|
|
|
||
|
|
if "`b2title'" == "" {
|
||
|
|
local b2title : variable label `theta'
|
||
|
|
if "`b2title'" == "" { local b2title "`theta'" }
|
||
|
|
}
|
||
|
|
if "`ci'" == "" { local b2title = "b2(s.e. of: `b2title')" }
|
||
|
|
else { local b2title = "b2(s.e. of: log[`b2title'])" }
|
||
|
|
|
||
|
|
if `"`t2title'"' == `""' {
|
||
|
|
local t2title = "Begg's funnel plot with pseudo"
|
||
|
|
local t2title = "`t2title' `level'% confidence limits"
|
||
|
|
local t2title "t2(`"`t2title'"')"
|
||
|
|
}
|
||
|
|
else if `"`t2title'"' == `"."' { local t2title }
|
||
|
|
else { local t2title "t2(`"`t2title'"')" }
|
||
|
|
|
||
|
|
if "`xlabel'" == "" { local xlabel "xla" }
|
||
|
|
else { local xlabel "xlabel(`xlabel')" }
|
||
|
|
|
||
|
|
if "`ylabel'" == "" { local ylabel "yla" }
|
||
|
|
else { local ylabel "ylabel(`ylabel')" }
|
||
|
|
|
||
|
|
if "`sort'" == "" { local sort "sort" }
|
||
|
|
|
||
|
|
if "`gap'" == "" { local gap "gap(3)" }
|
||
|
|
else { local gap "gap(`gap')" }
|
||
|
|
|
||
|
|
tempvar ll2 ul2 z mmm var w sw wl swl RRm
|
||
|
|
tempname oe
|
||
|
|
|
||
|
|
qui {
|
||
|
|
if `level' < 1 { local `level' =`level' * 100 }
|
||
|
|
local z = -1 * invnorm((1 - `level' / 100) / 2)
|
||
|
|
local obs1=_N+1
|
||
|
|
set obs `obs1'
|
||
|
|
replace `setheta'=0 in `obs1'
|
||
|
|
replace `theta' = . in `obs1'
|
||
|
|
gen `var' = `setheta'^2
|
||
|
|
gen `w' = 1/`var'
|
||
|
|
egen `sw' = sum(`w') if `touse'
|
||
|
|
gen `wl' = `w' * `theta'
|
||
|
|
egen `swl' = sum(`wl') if `touse'
|
||
|
|
gen `RRm' = `swl' / `sw'
|
||
|
|
scalar `oe' = `RRm'
|
||
|
|
egen `mmm' = min(`RRm')
|
||
|
|
replace `RRm' = `mmm' if `setheta' == 0
|
||
|
|
gen `ll2' = `RRm' - `z' * `setheta'
|
||
|
|
gen `ul2' = `RRm' + `z' * `setheta'
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`gweight'" != "" {
|
||
|
|
tempvar ww
|
||
|
|
qui gen `ww' = `w'
|
||
|
|
qui replace `ww' = 1 if `setheta' == 0
|
||
|
|
local gw "[w=`ww']"
|
||
|
|
}
|
||
|
|
else { local gw "" }
|
||
|
|
|
||
|
|
#delimit ;
|
||
|
|
graph `RRm' `ll2' `ul2' `theta' `setheta' `gw' if `touse',
|
||
|
|
`connect' `symbol' `t2title' `l2title' `b2title'
|
||
|
|
`l1title' `xlabel' `ylabel' `sort' `pen' `gap' `options';
|
||
|
|
#delimit cr
|
||
|
|
|
||
|
|
exit
|
||
|
|
end
|
||
|
|
|
||
|
|
* ***************************************************
|
||
|
|
|
||
|
|
program define egggph
|
||
|
|
version 6.0
|
||
|
|
|
||
|
|
* creates the Egger regression asymmetry plot graph
|
||
|
|
|
||
|
|
* Setup
|
||
|
|
|
||
|
|
syntax varlist(min=3 max=3) [if] [in] [, LEVel(integer 95) GWeight /*
|
||
|
|
*/ GAp(str) Connect(str) Symbol(str) SOrt Pen(str) T2title(str) /*
|
||
|
|
*/ B2title(str) YLAbel(str) XLAbel(str) L1title(str) L2title(str) *]
|
||
|
|
|
||
|
|
tokenize `varlist'
|
||
|
|
|
||
|
|
tempvar touse theta setheta
|
||
|
|
|
||
|
|
local theta `1'
|
||
|
|
local setheta `2'
|
||
|
|
local touse `3'
|
||
|
|
|
||
|
|
preserve
|
||
|
|
|
||
|
|
* Graph options
|
||
|
|
if "`connect'" == "" { local connect ".ll" }
|
||
|
|
else {
|
||
|
|
local lll = length("`connect'")
|
||
|
|
if `lll' == 1 { local connect "`connect'll" }
|
||
|
|
else if `lll' == 2 { local connect "`connect'l" }
|
||
|
|
}
|
||
|
|
local connect "co(`connect')"
|
||
|
|
|
||
|
|
if "`symbol'" == "" { local symbol "oid" }
|
||
|
|
else {
|
||
|
|
local lll = length("`symbol'")
|
||
|
|
if `lll' == 1 { local symbol "`symbol'id" }
|
||
|
|
else if `lll' == 2 { local symbol "`symbol'd" }
|
||
|
|
}
|
||
|
|
local symbol "sy(`symbol')"
|
||
|
|
|
||
|
|
if "`pen'" == "" { local pen "233" }
|
||
|
|
else {
|
||
|
|
local lll = length("`pen'")
|
||
|
|
if `lll' == 1 { local pen "`pen'33" }
|
||
|
|
else if `lll' == 2 {local pen = "`pen'" + substr("`pen'",2,1)}
|
||
|
|
}
|
||
|
|
local pen "pen(`pen')"
|
||
|
|
|
||
|
|
if "`l2title'" == "" { local l2title "standardized effect" }
|
||
|
|
local l2title "l2(`l2title')"
|
||
|
|
|
||
|
|
if "`l1title'" == "" { local l1title "" "" }
|
||
|
|
local l1title "l1(`l1title')"
|
||
|
|
|
||
|
|
if "`b2title'" == "" { local b2title "precision" }
|
||
|
|
local b2title = "b2(`b2title')"
|
||
|
|
|
||
|
|
if `"`t2title'"' == `""' {
|
||
|
|
local t2title = "Egger's publication bias plot"
|
||
|
|
local t2title "t2(`"`t2title'"')"
|
||
|
|
}
|
||
|
|
else if `"`t2title'"' == `"."' { local t2title t2(" ")}
|
||
|
|
else { local t2title "t2(`"`t2title'"')" }
|
||
|
|
|
||
|
|
|
||
|
|
if "`xlabel'" == "" { local xlabel "xla" }
|
||
|
|
else { local xlabel "xlabel(`xlabel')" }
|
||
|
|
|
||
|
|
if "`ylabel'" == "" { local ylabel "yla" }
|
||
|
|
else { local ylabel "ylabel(`ylabel')" }
|
||
|
|
|
||
|
|
if "`sort'" == "" { local sort "sort" }
|
||
|
|
|
||
|
|
if "`gap'" == "" { local gap "gap(3)" }
|
||
|
|
else { local gap "gap(`gap')" }
|
||
|
|
|
||
|
|
qui {
|
||
|
|
local obs1 = _N + 1
|
||
|
|
local obs2 = _N + 2
|
||
|
|
local obs3 = _N + 3
|
||
|
|
set obs `obs3'
|
||
|
|
replace `setheta' = 0 in `obs1'/`obs3'
|
||
|
|
replace `theta' = . in `obs1'/`obs3'
|
||
|
|
tempvar prec snd
|
||
|
|
gen `prec' = 1 / `setheta' if `setheta' > 0
|
||
|
|
gen `snd' = `theta' / `setheta' if `setheta' > 0
|
||
|
|
replace `prec' = 0 if `prec' == .
|
||
|
|
regr `snd' `prec' if `touse'
|
||
|
|
tempvar reg ci
|
||
|
|
capture matrix b = get(_b)
|
||
|
|
if _rc == 0 {
|
||
|
|
matrix V = get(VCE)
|
||
|
|
local df = e(N) - 2
|
||
|
|
gen `reg' = b[1,2] + `prec' * b[1,1]
|
||
|
|
gen `ci' = .
|
||
|
|
#delimit ;
|
||
|
|
replace `ci' = b[1,2] - sqrt(V[2,2])
|
||
|
|
* invt(`df', `level'/100) in `obs2' ;
|
||
|
|
replace `ci' = b[1,2] + sqrt(V[2,2])
|
||
|
|
* invt(`df', `level'/100) in `obs3' ;
|
||
|
|
#delimit cr
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
gen `reg' = .
|
||
|
|
gen `ci' = .
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`gweight'" != "" {
|
||
|
|
tempvar ww
|
||
|
|
qui gen `ww' = 1 / `setheta'^2
|
||
|
|
qui replace `ww' = 1 if `setheta' == 0
|
||
|
|
local gw "[w=`ww']"
|
||
|
|
}
|
||
|
|
else { local gw "" }
|
||
|
|
|
||
|
|
#delimit ;
|
||
|
|
graph `snd' `reg' `ci' `prec' `gw' if `touse', yli(0) xli(0)
|
||
|
|
`connect' `symbol' `t2title' `l2title' `b2title'
|
||
|
|
`l1title' `xlabel' `ylabel' `sort' `pen' `gap' `options';
|
||
|
|
#delimit cr
|
||
|
|
|
||
|
|
exit
|
||
|
|
end
|
||
|
|
|
||
|
|
* ***************************************************
|
||
|
|
|
||
|
|
*! ktau2 version 4.1.0 26sep97 TJS
|
||
|
|
program define ktau2
|
||
|
|
version 4.0
|
||
|
|
|
||
|
|
* modification of ktau to allow N==2, un-continuity-corrected
|
||
|
|
* z and p values, and to pass more parameters
|
||
|
|
|
||
|
|
local varlist "req ex min(2) max(2)"
|
||
|
|
local if "opt"
|
||
|
|
local in "opt"
|
||
|
|
parse "`*'"
|
||
|
|
parse "`varlist'", parse(" ")
|
||
|
|
local x "`1'"
|
||
|
|
local y "`2'"
|
||
|
|
tempname k N NN pval score se tau_a tau_b
|
||
|
|
tempname xt xt2 xt3 yt yt2 yt3
|
||
|
|
tempvar doit nobs order work
|
||
|
|
mark `doit' `in' `if'
|
||
|
|
markout `doit' `x' `y'
|
||
|
|
quietly count if `doit'
|
||
|
|
scalar `N' = _result(1)
|
||
|
|
if `N' < 2 { error 2001 }
|
||
|
|
local Nmac = `N'
|
||
|
|
qui {
|
||
|
|
gen long `order' = _n /* restore ordering at end */
|
||
|
|
replace `doit' = -`doit'
|
||
|
|
sort `doit' /* put obs for computation first */
|
||
|
|
gen double `work' = 0 /* using type double is fastest */
|
||
|
|
scalar `k' = 2
|
||
|
|
while (`k' <= `N') {
|
||
|
|
local kk = `k' - 1
|
||
|
|
#delimit ;
|
||
|
|
replace `work' = `work'
|
||
|
|
+ sign((`x' - `x'[`k'])*(`y' - `y'[`k']))
|
||
|
|
in 1/`kk' ; /* using "in" is fastest */
|
||
|
|
#delimit cr
|
||
|
|
scalar `k' = `k' + 1
|
||
|
|
}
|
||
|
|
replace `work' = sum(`work') in 1/`Nmac'
|
||
|
|
scalar `score' = `work'[`N']
|
||
|
|
/* Calculate ties on `x' */
|
||
|
|
egen long `nobs' = count(`x') in 1/`Nmac', by(`x')
|
||
|
|
tempvar nobsxm
|
||
|
|
egen `nobsxm' = max(`nobs')
|
||
|
|
/* Calculate correction term for ties on `x' */
|
||
|
|
replace `work' = sum((`nobs' - 1)*(2*`nobs' + 5)) in 1/`Nmac'
|
||
|
|
scalar `xt' = `work'[`N']
|
||
|
|
/* Calculate correction term for pairs of ties on `x' */
|
||
|
|
replace `work' = sum(`nobs' - 1) in 1/`Nmac'
|
||
|
|
scalar `xt2' = `work'[`N']
|
||
|
|
/* Calculate correction term for triplets of ties on `x' */
|
||
|
|
replace `work' = sum((`nobs' - 1)*(`nobs' - 2)) in 1/`Nmac'
|
||
|
|
scalar `xt3' = `work'[`N']
|
||
|
|
/* Calculate ties on `y' */
|
||
|
|
drop `nobs'
|
||
|
|
egen long `nobs' = count(`y') in 1/`Nmac', by(`y')
|
||
|
|
tempvar nobsym
|
||
|
|
egen `nobsym' = max(`nobs')
|
||
|
|
/* Calculate correction term for ties on `y' */
|
||
|
|
replace `work' = sum((`nobs' - 1)*(2*`nobs' + 5)) in 1/`Nmac'
|
||
|
|
scalar `yt' = `work'[`N']
|
||
|
|
/* Calculate correction term for pairs of ties on `y' */
|
||
|
|
replace `work' = sum(`nobs' - 1) in 1/`Nmac'
|
||
|
|
scalar `yt2' = `work'[`N']
|
||
|
|
/* Calculate correction term for triplets of ties on `y' */
|
||
|
|
replace `work' = sum((`nobs' - 1)*(`nobs' - 2)) in 1/`Nmac'
|
||
|
|
scalar `yt3' = `work'[`N']
|
||
|
|
/* Compute Kendall's tau-a, tau-b, s.e. of score, and pval */
|
||
|
|
scalar `NN' = `N'*(`N' - 1)
|
||
|
|
scalar `tau_a' = 2*`score'/`NN'
|
||
|
|
scalar `tau_b' = 2*`score'/sqrt((`NN' - `xt2')*(`NN' - `yt2'))
|
||
|
|
#delimit ;
|
||
|
|
scalar `se' = `NN'*(2*`N' + 5);
|
||
|
|
tempname tmax;
|
||
|
|
scalar `tmax' = max(`nobsxm', `nobsym');
|
||
|
|
if `tmax' > 1 { scalar `se' = `se'
|
||
|
|
- (`xt' - `yt')
|
||
|
|
+ `xt3'*`yt3'/(9*`NN'*(`N' - 2))
|
||
|
|
+ `xt2'*`yt2'/(2*`NN') } ;
|
||
|
|
scalar `se' = sqrt((1/18)*`se');
|
||
|
|
#delimit cr
|
||
|
|
local zcc = (abs(`score') - 1) / `se'
|
||
|
|
local z = `score '/ `se'
|
||
|
|
tempname pvalcc
|
||
|
|
if `score' == 0 {
|
||
|
|
scalar `pval' = 1
|
||
|
|
scalar `pvalcc' = 1
|
||
|
|
}
|
||
|
|
else scalar `pvalcc' = 2*(1 - normprob((abs(`score') - 1)/`se'))
|
||
|
|
else scalar `pval' = 2*(1 - normprob(abs(`score')/`se'))
|
||
|
|
/* Restore original ordering of data set */
|
||
|
|
sort `order'
|
||
|
|
}
|
||
|
|
/* Print results */
|
||
|
|
#delimit ;
|
||
|
|
di _n
|
||
|
|
in gr " Number of obs = " in ye %7.0f `N' _n
|
||
|
|
in gr "Kendall's tau-a = " in ye %12.4f `tau_a' _n
|
||
|
|
in gr "Kendall's tau-b = " in ye %12.4f `tau_b' _n
|
||
|
|
in gr "Kendall's score = " in ye %7.0f `score' _n
|
||
|
|
in gr " SE of score = " in ye %11.3f `se' _c ;
|
||
|
|
if `xt2' > 0 | `yt2' > 0 { di in gr " (corrected for ties)" _c } ;
|
||
|
|
di _n(2)
|
||
|
|
in gr "Test of Ho: `x' and `y' independent" _n
|
||
|
|
in gr " z = " in ye %12.2f `z' _n
|
||
|
|
in gr " Pr > |z| = " in ye %12.4f = `pval' _n(2)
|
||
|
|
in gr " z = " in ye %12.2f sign(`score')*`zcc' _n
|
||
|
|
in gr " Pr > |z| = " in ye %12.4f = `pvalcc'
|
||
|
|
in gr " (continuity corrected)" ;
|
||
|
|
#delimit cr
|
||
|
|
local c = 0
|
||
|
|
if `xt2' > 0 | `yt2' > 0 { local c = 1 }
|
||
|
|
global S_1 = `N'
|
||
|
|
global S_2 = `tau_a'
|
||
|
|
global S_3 = `tau_b'
|
||
|
|
global S_4 = `score'
|
||
|
|
global S_5 = `se'
|
||
|
|
global S_6 = `pval'
|
||
|
|
global S_7 = `z'
|
||
|
|
global S_8 = `pvalcc'
|
||
|
|
global S_9 = `zcc'
|
||
|
|
global S_10 = `c'
|
||
|
|
end
|
||
|
|
|
||
|
|
* ***************************************************
|
||
|
|
|
||
|
|
*! ifexp version 1.2.1 19nov98 TJS
|
||
|
|
program define ifexp
|
||
|
|
version 5.0
|
||
|
|
|
||
|
|
while substr("`1'",1,2) != "if" {
|
||
|
|
macro shift
|
||
|
|
}
|
||
|
|
|
||
|
|
local if "required"
|
||
|
|
local options "noVARlabel noVALlabel noUNabbrev noSPace Color(str)"
|
||
|
|
parse "`*'"
|
||
|
|
parse "`if'", parse(" ><=&|!()")
|
||
|
|
|
||
|
|
if "`color'" == "" {local color = "b"} /* set display color */
|
||
|
|
local color = lower(substr("`color'",1,1))
|
||
|
|
if index("wbgyr","`color'") == 0 {local color = "b"}
|
||
|
|
|
||
|
|
if "`if'" != "" {
|
||
|
|
local ifst " if " /* handle leading "if " */
|
||
|
|
local lif = 5
|
||
|
|
macro shift
|
||
|
|
|
||
|
|
if "`space'" != "nospace" { /* set up space */
|
||
|
|
local sp " "
|
||
|
|
local sp1 = 1
|
||
|
|
}
|
||
|
|
else { local sp1 = 0 }
|
||
|
|
|
||
|
|
di
|
||
|
|
while "`1'" != "" { /* start main loop on tokens */
|
||
|
|
capture confirm variable `1' /* var name ? */
|
||
|
|
local rc _rc
|
||
|
|
|
||
|
|
local m = 0
|
||
|
|
local n = length("`1'")
|
||
|
|
if `n' + `lif' > 78 { /* no room for token ? */
|
||
|
|
di in `color' "`ifst'"
|
||
|
|
local ifst "> "
|
||
|
|
local lif = 5
|
||
|
|
local n = 0
|
||
|
|
}
|
||
|
|
/* is punctuation */
|
||
|
|
if "`1'" == "==" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == "<=" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == ">=" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == "!=" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == "<" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == ">" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == "=" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == "&" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == "|" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == "!" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == "(" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
else if "`1'" == ")" { local ifst "`ifst'`1'`sp'" }
|
||
|
|
|
||
|
|
else if `rc' == 0 { /* is a variable name */
|
||
|
|
local tif "`1'"
|
||
|
|
if "`unabbrev'" != "nounabbrev" { /* allow unabbrev */
|
||
|
|
unabbrev `1'
|
||
|
|
local tif "$S_1"
|
||
|
|
}
|
||
|
|
if "`varlabel'" != "novarlabel" { /* allow var label ? */
|
||
|
|
local vrl : variable label `1'
|
||
|
|
if "`vrl'" != "" { local tif "`vrl'" } /* var label ? */
|
||
|
|
}
|
||
|
|
|
||
|
|
local wc: word count `tif' /* process > 1 word in label */
|
||
|
|
local w 1
|
||
|
|
while `w' < `wc' {
|
||
|
|
local wif: word `w' of `tif'
|
||
|
|
local wm = length("`wif'")
|
||
|
|
if `wm' + `lif' > 78 {
|
||
|
|
di in `color' "`ifst'"
|
||
|
|
local ifst "> "
|
||
|
|
local lif = 5
|
||
|
|
local n = `m'
|
||
|
|
local m = 0
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
local lif = `lif' + `wm' + `sp1'
|
||
|
|
local ifst "`ifst'`wif'`sp'"
|
||
|
|
local w = `w' + 1
|
||
|
|
}
|
||
|
|
}
|
||
|
|
local tif: word `wc' of `tif' /* do last or only word */
|
||
|
|
|
||
|
|
local m = length("`tif'") /* no room for token ? */
|
||
|
|
if `m' + `lif' > 78 {
|
||
|
|
di in `color' "`ifst'"
|
||
|
|
local ifst "> "
|
||
|
|
local lif = 5
|
||
|
|
local n = `m'
|
||
|
|
local m = 0
|
||
|
|
}
|
||
|
|
local ifst "`ifst'`tif'`sp'"
|
||
|
|
local vll : value label `1'
|
||
|
|
}
|
||
|
|
|
||
|
|
else { /* is a value */
|
||
|
|
local tif "`1'" /* allow value label ? */
|
||
|
|
if "`vll'" != "" & "`vallabel'" != "novallabel" {
|
||
|
|
local lv : label `vll' `1'
|
||
|
|
if "`lv'" != "" { local tif "`lv'" } /* value label ? */
|
||
|
|
}
|
||
|
|
|
||
|
|
local wc: word count `tif' /* process > 1 word in label */
|
||
|
|
local w 1
|
||
|
|
while `w' < `wc' {
|
||
|
|
local wif: word `w' of `tif'
|
||
|
|
local wm = length("`wif'")
|
||
|
|
if `wm' + `lif' > 78 {
|
||
|
|
di in `color' "`ifst'"
|
||
|
|
local ifst "> "
|
||
|
|
local lif = 5
|
||
|
|
local n = `m'
|
||
|
|
local m = 0
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
local lif = `lif' + `wm' + `sp1'
|
||
|
|
local ifst "`ifst'`wif'`sp'"
|
||
|
|
local w = `w' + 1
|
||
|
|
}
|
||
|
|
}
|
||
|
|
local tif: word `wc' of `tif' /* do last or only word */
|
||
|
|
|
||
|
|
local m = length("`tif'") /* no room for token ? */
|
||
|
|
if `m' + `lif' > 78 {
|
||
|
|
di in `color' "`ifst'"
|
||
|
|
local ifst "> "
|
||
|
|
local lif = 5
|
||
|
|
local n = `m'
|
||
|
|
local m = 0
|
||
|
|
}
|
||
|
|
local ifst "`ifst'`tif'`sp'"
|
||
|
|
}
|
||
|
|
/* set new line length */
|
||
|
|
if `m' > 0 { local lif = `lif' + `m' + `sp1' }
|
||
|
|
else { local lif = `lif' + `n' + `sp1' }
|
||
|
|
macro shift /* next token */
|
||
|
|
} /* end main loop on tokens */
|
||
|
|
|
||
|
|
di in `color' "`ifst'" /* display expanded if */
|
||
|
|
|
||
|
|
global S_1 "`ifst'"
|
||
|
|
exit
|
||
|
|
end
|