Bookmark and Share

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/


© Copyright 1996–2018 StataCorp LLC   |   Terms of use   |   Privacy   |   Contact us   |   Site index