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: "Can Your Results be Replicated?" (Stata error?)
From
Joerg Luedicke <[email protected]>
To
[email protected]
Subject
Re: st: "Can Your Results be Replicated?" (Stata error?)
Date
Fri, 13 Sep 2013 11:55:09 -0400
I wonder how this became a "Stata vs. R" thing? Perhaps the
"R-blogger" mentioned by Phil in his original posting simply assumed
that the problem would not have occurred if the authors of the
original 2009 study would have used R (that is, some R package)
instead of Stata. However, using the same data as in my Stata example
above, the -geeglm- function from the -geepack- package (version
1.1-6) also provides coefficient estimates etc.:
> #--------------------------------
> require(Hmisc)
> require(geepack)
>
> dat = stata.get("gee_check1.dta")
>
> M1 <- geeglm( y ~ x.1 + x.3, data=dat, id=id,
+ family=binomial(link="logit"),
+ corstr="exchangeable")
> summary(M1)
Call:
geeglm(formula = y ~ x.1 + x.3, family = binomial(link = "logit"),
data = dat, id = id, corstr = "exchangeable")
Coefficients:
Estimate Std.err Wald Pr(>|W|)
(Intercept) 4.75e-01 1.18e-01 1.63e+01 5.4e-05 ***
x.1 -1.73e+07 3.80e+04 2.08e+05 < 2e-16 ***
x.3 -2.05e-01 1.50e-01 1.86e+00 0.17
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Estimated Scale Parameters:
Estimate Std.err
(Intercept) 0.667 0.0134
Correlation: Structure = exchangeable Link = identity
Estimated Correlation Parameters:
Estimate Std.err
alpha 0.0187 0.0161
Number of clusters: 100 Maximum cluster size: 10
> #--------------------------------
If anything, this underlines the importance that users make reasonable
modeling choices as one can surely fit a lot of garbage models without
the computer program at hand complaining about it.
Joerg
On Fri, Sep 13, 2013 at 10:54 AM, Stas Kolenikov <[email protected]> wrote:
> Anders, thanks for the rigorous follow up.
>
> -- Stas Kolenikov, PhD, PStat (ASA, SSC)
> -- Senior Survey Statistician, Abt SRBI
> -- Opinions stated in this email are mine only, and do not reflect the
> position of my employer
> -- http://stas.kolenikov.name
>
>
>
> On Fri, Sep 13, 2013 at 9:51 AM, Anders Alexandersson
> <[email protected]> wrote:
>> Joerg is correct, and I agree with Stas. I contacted the author Mark
>> Bell who replied:
>>
>> "We never claim the error was in the firth logit command. The problem
>> was that Rauchhaus was using a logit GEE as a way to deal with
>> separation (xtgee). GEE is not a way to deal with separation, and the
>> model should have been unidentified, but stata nonetheless returned a
>> coefficient on the variable causing separation. This is the error we
>> were referring to."
>>
>> Anders
>>
>> On Fri, Sep 13, 2013 at 10:15 AM, Stas Kolenikov <[email protected]> wrote:
>>> Maybe it's time Stata Corp picks up -firthlogit-, solidifies it and
>>> makes it an official command.
>>>
>>> -- Stas Kolenikov, PhD, PStat (ASA, SSC)
>>> -- Senior Survey Statistician, Abt SRBI
>>> -- Opinions stated in this email are mine only, and do not reflect the
>>> position of my employer
>>> -- http://stas.kolenikov.name
>>>
>>>
>>>
>>> On Fri, Sep 13, 2013 at 9:01 AM, Joerg Luedicke
>>> <[email protected]> wrote:
>>>> After having a quick glance at their paper
>>>> (http://jcr.sagepub.com/content/early/2013/08/19/0022002713499718.abstract?papetoc)
>>>> it seems that they are talking about a problem with Stata's -xtgee-
>>>> command which, in the case of separation in a logit model, provides
>>>> nonsense results as opposed to omitting predictors or the like. Below
>>>> is a toy example showing what seems to be the problem. However,
>>>> finding an effect of something like "x is 3 million times less likely
>>>> than y" and not getting suspicious rather looks like sloppy research
>>>> to me in the first place.
>>>>
>>>> Joerg
>>>>
>>>>
>>>> *------------------------
>>>> clear
>>>> set obs 100
>>>> set seed 123
>>>>
>>>> gen id = _n
>>>> gen ui = rnormal(0,0.5)
>>>>
>>>> expand 10
>>>> bys id : gen year = _n
>>>> gen x = cond(mod(_n-1, 3) == 1, 1, cond(mod(_n-1, 3) == 0, 2, 3))
>>>> tab x, gen(x_)
>>>>
>>>> gen xb = 1 / (1 + exp(-(0.3*x_2 + 0.3*x_3 + ui)))
>>>> gen y = rbinomial(1,xb)
>>>> replace y = 0 if x_1 == 1
>>>> tab y x
>>>>
>>>> xtset id year
>>>> xtgee y i.x, fam(binomial) link(logit)
>>>> melogit y i.x || id:
>>>> logit y i.x
>>>> *------------------------
>>>>
>>>> On Fri, Sep 13, 2013 at 10:41 AM, Richard Williams
>>>> <[email protected]> wrote:
>>>>> At 08:14 AM 9/13/2013, Anders Alexandersson wrote:
>>>>>>
>>>>>> I just made this reply on the blog:
>>>>>>
>>>>>> "Where is the error in Stata? The author’s so called “Do-File for
>>>>>> Analyses.txt” is actually not a Stata do file but it does refer to
>>>>>> Stata’s user-written command -firthlogit- from SSC. Please provide a
>>>>>> reproducible do-file in Stata.The claim that results and conclusions
>>>>>> were due to an error in Stata is not supported."
>>>>>> See
>>>>>> http://politicalsciencereplication.wordpress.com/2013/09/11/guest-blog-how-to-persuade-journals-to-accept-your-replication-paper/comment-page-1/#comment-653
>>>>>
>>>>>
>>>>> Even if -firthlogit- did get it wrong, it is a bit of a stretch to imply
>>>>> that Stata has some terrible flaw. Stata Corp can hardly be held responsible
>>>>> for flaws in programs it did not write.
>>>>>
>>>>> When I first sent a program to SSC, I thought there might be some sort of
>>>>> exhaustive review process before it was released to the public. I got the
>>>>> feeling that wasn't the case when I got a message less than an hour later
>>>>> saying the program had been posted. Most user-written routines are fine but
>>>>> people should realize they haven't undergone the kind of testing that
>>>>> official programs have. And even in this case, we don't have any proof yet
>>>>> that firthlogit did get it wrong.
>>>>>
>>>>>
>>>>>> Anders Alexandersson
>>>>>> [email protected]
>>>>>>
>>>>>> On Fri, Sep 13, 2013 at 8:10 AM, Philip Jones
>>>>>> <[email protected]> wrote:
>>>>>> > Hi all,
>>>>>> >
>>>>>> > I found a link on my Twitter feed this AM, purporting to show how
>>>>>> > Stata "made a mistake" that R did not make:
>>>>>> >
>>>>>> > http://www.r-bloggers.com/can-your-results-be-replicated/
>>>>>> >
>>>>>> > which actually points to:
>>>>>> >
>>>>>> >
>>>>>> > http://politicalsciencereplication.wordpress.com/2013/09/11/guest-blog-how-to-persuade-journals-to-accept-your-replication-paper/
>>>>>> >
>>>>>> > I realize that "r-bloggers" is likely not the most bias-free site when
>>>>>> > it comes to reviewing/rating stats packages, but has anyone got an
>>>>>> > idea as to what is actually going on here? Is Stata really at fault?
>>>>>> >
>>>>>> > Regards,
>>>>>> >
>>>>>> > Phil
>>>>>> > @pmgjones on Twitter
>>>>>> > *
>>>>>> > * For searches and help try:
>>>>>> > * http://www.stata.com/help.cgi?search
>>>>>> > * http://www.stata.com/support/faqs/resources/statalist-faq/
>>>>>> > * http://www.ats.ucla.edu/stat/stata/
>>>>>>
>>>>>> *
>>>>>> * For searches and help try:
>>>>>> * http://www.stata.com/help.cgi?search
>>>>>> * http://www.stata.com/support/faqs/resources/statalist-faq/
>>>>>> * http://www.ats.ucla.edu/stat/stata/
>>>>>
>>>>>
>>>>> -------------------------------------------
>>>>> Richard Williams, Notre Dame Dept of Sociology
>>>>> OFFICE: (574)631-6668, (574)631-6463
>>>>> HOME: (574)289-5227
>>>>> EMAIL: [email protected]
>>>>> WWW: http://www.nd.edu/~rwilliam
>>>>>
>>>>>
>>>>>
>>>>> *
>>>>> * For searches and help try:
>>>>> * http://www.stata.com/help.cgi?search
>>>>> * http://www.stata.com/support/faqs/resources/statalist-faq/
>>>>> * http://www.ats.ucla.edu/stat/stata/
>>>>
>>>> *
>>>> * For searches and help try:
>>>> * http://www.stata.com/help.cgi?search
>>>> * http://www.stata.com/support/faqs/resources/statalist-faq/
>>>> * http://www.ats.ucla.edu/stat/stata/
>>>
>>> *
>>> * For searches and help try:
>>> * http://www.stata.com/help.cgi?search
>>> * http://www.stata.com/support/faqs/resources/statalist-faq/
>>> * http://www.ats.ucla.edu/stat/stata/
>>
>> *
>> * For searches and help try:
>> * http://www.stata.com/help.cgi?search
>> * http://www.stata.com/support/faqs/resources/statalist-faq/
>> * http://www.ats.ucla.edu/stat/stata/
>
> *
> * For searches and help try:
> * http://www.stata.com/help.cgi?search
> * http://www.stata.com/support/faqs/resources/statalist-faq/
> * http://www.ats.ucla.edu/stat/stata/
*
* For searches and help try:
* http://www.stata.com/help.cgi?search
* http://www.stata.com/support/faqs/resources/statalist-faq/
* http://www.ats.ucla.edu/stat/stata/