Berechnung des Schwerkraftmodells in R- und Stata-Software: Warum sind Koeffizienten gleich, aber Standardfehler unterschiedlich?

Wir haben das Schwerkraftmodell in der R- und Stata-Software berechnet.

Für Berechnungen haben wir das Standardpaket glmmin R (mit Parameter family = quasipoisson) und ppmlin Stata verwendet.

Rufen Sie das Berechnungsverfahren in R auf:

summary(glmm<-glm(formula=exports ~ ln_GDPimporter + ln_GDPexporter + 
    ln_GDPimppc + ln_GDPexppc + ln_Distance + ln_Tariff + ln_ExchangeRate +
    Contig + Comlang + Colony_CIS + EAEU_CIS + EU_European_Union, 
    family=quasipoisson(link="log"),data=data_pua))

Die Ergebnisse in R waren:

Mit denselben Daten haben wir in Stata folgende Berechnungen ppmldurchgeführt:

ppml exports ln_gdpimporter ln_gdpexporter ln_gdpimppc ln_gdpexppc ln_distance ln_tariff ln_exchangerate contig comlang colony_cis eaeu_cis eu_european_union

Die Ergebnisse der Berechnungen in Stata waren wie folgt:

Wie Sie sehen können, sind die Modellkoeffizienten (zweite Spalte in der Ergebnistabelle) mindestens bis zur Dezimalstelle der 4. Markierung gleich.

Andere Ergebnisse (aus der dritten Spalte in der Ergebnistabelle) sind jedoch nicht dieselben.

• Können Sie Unterschiede in den Ergebnissen erklären?

• Warum sind die Koeffizienten insbesondere gleich (die ersten Spalten der Ergebnistabelle), Standardfehler jedoch nicht?

Ich stelle fest, dass in der Stata-Koeffiziententabelle "Robust Std. Err." Erwähnt wurde, während glmmwahrscheinlich keine robusten Fehler verwendet werden. Das würde SE-Unterschiede erklären.

Auch ppmlscheint tatsächlich „nicht signifikant“ Regressoren fallen, und R quasipoissonFamilie ermöglicht über Dispersion in einer Weise , die anders ist, sagen wir, negative binomiale Regression, die vielleicht anders ppml.

Ich habe festgestellt, dass Sie an einigen Stellen gefragt haben, welches R-Paket ppmlfür (ökonomische) Schwerkraftmodelle gleichwertige Ergebnisse liefern würde , und keine Antworten erhalten haben. Es tut mir leid, das zu sehen und ich wünschte, ich könnte eine fundiertere Empfehlung geben. Es scheint, dass Sie eine Poisson-Regression mit robusten Standardfehlern benötigen, die Nullwerte verarbeitet. Ich bin mir nicht sicher, welche R-Pakete das unterstützen. (Nicht sicher, ob ppmlGriffe über Dispersion oder nicht.)

Bayesianische Regressionspakete wie rstanarmHeteroskedastizität könnten robuster sein, aber ich bin mir nicht sicher. Ich würde eher so etwas wie eine student_tFamilie dafür benutzen , aber du musst es benutzen, poissondamit ich mir der Antwort dort nicht sicher bin. Sie können die negative Binomialfamilie ( neg_binomial_2in rstanarm's stan_glm) ausprobieren , die auch die Überdispersion behandelt und möglicherweise robuster ist als quasipoisson.

Siehe auch: Wann sollten robuste Standardfehler in der Poisson-Regression verwendet werden?

Um Waynes ausgezeichnete Antwort zu erweitern, ppmlwird eine robuste (bis heteroskedastische) Varianz-Kovarianz-Matrix sowie eine endliche Stichprobenanpassung an diese Matrix verwendet, um die Verzerrung zu verringern.

Diese sind sehr ähnlich zu dem, was sandwich()aus dem gleichnamigen Paket in R berechnet wird. Der einzige Unterschied besteht darin, wie die Anpassung der endlichen Stichprobe durchgeführt wird. In der sandwich(...)Funktion wird standardmäßig überhaupt keine Anpassung der endlichen Stichprobe vorgenommen, dh das Sandwich wird durch 1 / n geteilt, wobei n die Anzahl der Beobachtungen ist. Alternativ sandwich(..., adjust = TRUE)kann verwendet werden, was durch 1 / (n - k) dividiert, wobei k die Anzahl der Regressoren ist. Stata teilt jedoch durch 1 / (n - 1).

So können Sie R dazu bringen, mit Stata übereinzustimmen, indem Sie eine benutzerdefinierte Sandwich-Varianz mit einem Anpassungsfaktor von 1 / (n-1) verwenden:

. clear

. set more off

. capture ssc install rsource

. use http://personal.lse.ac.uk/tenreyro/mock, clear

. saveold ~/Desktop/mock, version(12) replace
(saving in Stata 12 format, which can be read by Stata 11 or 12)
file ~/Desktop/mock.dta saved

. rsource, terminator(XXX) rpath("/usr/local/bin/R") roptions("--vanilla")
Assumed R program path: "/usr/local/bin/R"

Loading required package: zoo

Attaching package: 'zoo'

The following objects are masked from 'package:base':

    as.Date, as.Date.numeric

Beginning of R output

R version 3.2.4 (2016-03-10) -- "Very Secure Dishes"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-apple-darwin13.4.0 (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

  Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

>   library("foreign")
>   library("sandwich")
>   library("lmtest")
>   mock<-read.dta("~/Desktop/mock.dta")
>   glmm<-glm(formula=y ~ x + w, family=quasipoisson(link="log"),data=mock)
> 
>   sandwich1 <- function(object, ...) sandwich(object) * nobs(object) / (nobs(object) - 1)
>   coeftest(glmm,vcov=sandwich1)

z test of coefficients:

            Estimate Std. Error z value  Pr(>|z|)    
(Intercept) 0.516969   0.098062  5.2718 1.351e-07 ***
x           0.125657   0.101591  1.2369    0.2161    
w           0.013410   0.710752  0.0189    0.9849    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

> 
End of R output

. 
. ppml y x w

note: checking the existence of the estimates

Number of regressors excluded to ensure that the estimates exist: 0
Number of observations excluded: 0

note: starting ppml estimation
note: y has noninteger values

Iteration 1:   deviance =  139.7855
Iteration 2:   deviance =  137.7284
Iteration 3:   deviance =  137.7222
Iteration 4:   deviance =  137.7222

Number of parameters: 3
Number of observations: 100
Pseudo log-likelihood: -173.89764
R-squared: .01628639
Option strict is: off
------------------------------------------------------------------------------
             |               Robust
           y |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
           x |   .1256565   .1015913     1.24   0.216    -.0734588    .3247718
           w |   .0134101   .7107518     0.02   0.985    -1.379638    1.406458
       _cons |   .5169689   .0980624     5.27   0.000     .3247702    .7091676
------------------------------------------------------------------------------

Hier ist der Stata / R-Code, der die obige Ausgabe generiert. Ich verwende rsource, um R von Stata aus auszuführen (und Sie müssen das rpath()Folgende anpassen , um es an Ihr Setup anzupassen), aber das ist nicht wirklich notwendig: Sie können das rsourceTeil einfach von R aus ausführen .

clear
set more off
capture ssc install rsource

use http://personal.lse.ac.uk/tenreyro/mock, clear
saveold ~/Desktop/mock, version(12) replace

rsource, terminator(XXX) rpath("/usr/local/bin/R") roptions("--vanilla")
  library("foreign")
  library("sandwich")
  library("lmtest")
  mock<-read.dta("~/Desktop/mock.dta")
  glmm<-glm(formula=y ~ x + w, family=quasipoisson(link="log"),data=mock)

  sandwich1 <- function(object, ...) sandwich(object) * nobs(object) / (nobs(object) - 1)
  coeftest(glmm,vcov=sandwich1)  
XXX 

ppml y x w