Notice: On April 23, 2014, Statalist moved from an email list to a forum, based at statalist.org.
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: st: RE: margeff after the glm command
From
Sabrina Carrossa <[email protected]>
To
[email protected]
Subject
Re: st: RE: margeff after the glm command
Date
Thu, 13 May 2010 17:24:46 +0100
Ok Martin,
I thought that the problem was about executing the margeff command
after the glm one, but obviously this is not the case.
I just tried to execute the margeff command after the "firthlogit
ptogram" (you can find the do file here:
http://www.stata.com/statalist/archive/2005-11/msg00786.html), but
stata shows " __000006 not found r(111) ".
To be more clear, I'd better paste all the commands:
***************************************************************
* firstly I execute the following program written by Joseph Coveney
***************************************************************
clear
set more off
input byte uti byte GrpSize byte age byte oc byte vic byte ///
vicl byte vis byte dia
1 2 0 0 0 0 1 0
14 16 0 0 1 0 0 0
3 4 0 0 1 0 1 0
10 18 0 0 1 1 0 0
12 30 0 0 1 1 1 0
1 1 0 0 1 1 1 1
44 86 0 1 0 0 0 0
3 3 0 1 0 0 0 1
0 1 0 1 0 0 1 0
15 16 0 1 1 0 0 0
1 1 0 1 1 0 0 1
2 2 0 1 1 0 1 0
7 12 0 1 1 1 0 0
3 9 0 1 1 1 1 0
2 2 1 0 1 0 0 0
1 1 1 0 1 0 1 0
1 1 1 0 1 0 1 1
1 3 1 0 1 1 0 0
0 4 1 0 1 1 1 0
1 1 1 0 1 1 1 1
5 19 1 1 0 0 0 0
0 1 1 1 0 0 1 0
3 4 1 1 1 0 0 0
0 2 1 1 1 1 1 0
end
rename uti count1
generate byte count0 = GrpSize - count1
drop GrpSize
reshape long count, i(age-dia) j(positive)
drop if count == 0
expand count
drop count
rename oc oral
rename vic condom
rename vicl lubri
rename vis sperm
rename dia diaphrag
replace age = !age // Heinze & Ploner flipped this
*
capture program drop firthlogit
program define firthlogit, eclass sortpreserve byable(recall)
syntax varlist(min = 1) [if] [in], [or]
tempvar hat key split false pseudo_response weight
preserve
quietly {
glm `varlist' , family(binomial) link(logit) `nolog'
local depvar `e(depvar)'
local indepvarlist = ///
subinstr("`varlist'", "`depvar'", "", 1)
predict `hat', hat
* Splitting
generate long `key' = _n
generate byte `split' = 2
expand `split'
drop `split'
bysort `key': generate byte `false' = _n - 1
generate byte `pseudo_response' = `depvar'
bysort `key' (`false'): replace `pseudo_response' = ///
1 - `depvar' if _n == _N
generate float `weight' = 1 + `hat' / 2
replace `weight' = `hat' / 2 if `false'
}
* Loop to convergence
local last_ll `e(ll)'
local criterion = 1
local i = 1
while `criterion' > 10^-4 {
quietly {
glm `pseudo_response' `indepvarlist' ///
[iweight = `weight'], family(binomial) link(logit)
local criterion = abs((`last_ll' - e(ll) ) / `last_ll')
drop `hat'
predict `hat', hat
replace `weight' = 1 + `hat' / 2
replace `weight' = `hat' / 2 if `false'
}
display in smcl as text "Firth", "Last Log-", ///
"This Log-", "Fractional"
display in smcl as text "Iteration", "likelihood", ///
"likelihood", "Change"
display in smcl as result "`++i'", `last_ll', ///
`e(ll)', `criterion'
local last_ll = e(ll)
}
if "`or'" != "" {
glm , eform
}
else {
glm
}
end
*
firthlogit positive age-diaphrag
exit
***************************************************************
* then I try to execute the margeff command
***************************************************************
margeff, dummies( age oral condom lubri sperm diaphrag) replace
/* but I got the error message: "__000006 not found
r(111); */
Hope this is more clear,
sabry
*
* 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/