[R] About calculation of the gravity model in R and STATA software

Jeff Newmiller jdnewmil at dcn.davis.ca.us
Tue Mar 15 20:37:57 CET 2016


I am not the person to answer your question, but have some suggestions:

Make your examples reproducible so others can confirm your results and explore other issues you may not have seen. [1] [2] 

Post your question on the R-sig-mixed mailing list, where mixed models experts hang out, instead of here where the R language is on topic rather than contributed packages like glmm.

I don't know SAS or the glmm package, but the fact that the top description line for the glmm package says it uses Monte  Carlo suggests to me that it might not be "standard" in SAS or R and that your results may vary depending on how you configure it and how the random number seed is initialized. 

[1] http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example
[2] See also the Posting Guide mentioned at the bottom of this and every email on this list. 
-- 
Sent from my phone. Please excuse my brevity.

On March 15, 2016 12:02:07 PM PDT, "Сергей С." <salnsg at gmail.com> wrote:
>Dear colleagues!
>
>We spent calculation of the gravity model in R and STATA software.
>For calculations we used the standard package glmm in R (with parameter
>family = quasipoisson)
>and ppml in STATA.
>
>Call the calculation procedure in R:
>
>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))
>
>The results of the calculations in R following:
>
>------------------------------------------------------------------------
>
>Coefficients:
>                           Estimate    Std. Error  t value Pr(>|t|)
>(Intercept)           -12.53224   15.30072  -0.819  0.41357
>ln_GDPimporter       0.10180    0.14988   0.679   0.49765
>ln_GDPexporter       0.14612    0.79823    0.183   0.85491
>ln_GDPimppc           0.34998    0.30247   1.157   0.24840
>ln_GDPexppc           0.65811    0.82189    0.801   0.42409
>ln_Distance             0.21838    0.16623    1.314   0.19020
>ln_Tariff                 -0.05499    0.04913  -1.119  0.26411
>ln_ExchangeRate     -0.11748    0.04275  -2.748  0.00646 **
>Contig                      1.48321    0.28684   5.171  4.92e-07 ***
>Comlang                  1.50727    0.26199    5.753   2.67e-08 ***
>Colony_CIS               2.15272    0.46899   4.590  7.16e-06 ***
>EAEU_CIS                -0.94417    0.29315  -3.221  0.00146 **
>EU_European_Union  -0.08335    0.76733  -0.109  0.91359
>---
>Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
>(Dispersion parameter for quasipoisson family taken to be 2100.979)
>
>    Null deviance: 1886758  on 251  degrees of freedom
>Residual deviance:  316332  on 239  degrees of freedom
>AIC: NA
>
>Number of Fisher Scoring iterations: 8
>
>------------------------------------------------------------------------
>
>On the same data, we done calculations in STATA, using ppml procedure.
>Call of the calculation procedure in STATA was next:
>
>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
>
>The results of the calculations in STATA were following:
>
>------------------------------------------------------------------------
>
>Iteration 1:   deviance =  425911.3
>Iteration 2:   deviance =  327020.8
>Iteration 3:   deviance =  316763.3
>Iteration 4:   deviance =  316335.1
>Iteration 5:   deviance =  316332.3
>Iteration 6:   deviance =  316332.3
>Iteration 7:   deviance =  316332.3
>
>Number of parameters: 13
>Number of observations: 252
>Pseudo log-likelihood: -158930.64
>R-squared: .75348104
>Option strict is: off
>-----------------------------------------------------------------------------------
>                           |                    Robust
>          exports       |      Coef.      Std. Err.          z    P>|z|
>   [95% Conf. Interval]
>------------------
>+----------------------------------------------------------------
>   ln_gdpimporter    |   .1018021   .0982091     1.04   0.300
>-.0906843    .2942885
>   ln_gdpexporter    |   .1461135   1.084255     0.13   0.893
>-1.978988    2.271215
>      ln_gdpimppc     |    .349982    .201011      1.74   0.082
>-.0439924    .7439564
>      ln_gdpexppc     |   .6581201   1.098236     0.60   0.549
>-1.494383    2.810624
>      ln_distance       |   .2183809    .156757     1.39    0.164
>-.0888572     .525619
>        ln_tariff          |  -.0549914   .0551489    -1.00   0.319
>-.1630811    .0530984
>  ln_exchangerate    |  -.1174816   .0343881    -3.42   0.001
>-.1848812   -.0500821
>           contig          |   1.483213    .168467     8.80   0.000
>1.153024    1.813402
>          comlang       |   1.507272   .2745761     5.49   0.000
>.9691126    2.045431
>       colony_cis        |   2.152723   .2338133     9.21   0.000
>1.694457    2.610988
>         eaeu_cis        |  -.9441651   .2469764    -3.82   0.000
>-1.42823   -.4601003
>eu_european_union |  -.0833477   .4955678    -0.17   0.866    
>-1.054643
>.8879474
>            _cons         |  -12.53206   21.18599    -0.59   0.554
>-54.05585    28.99172
>-----------------------------------------------------------------------------------
>
>As you can see, model coefficients (second column in the results table)
>are
>the same at least until the 4th mark (!)
>However, other results (columns in the table of results, since the
>third)
>is not the same.
>Could you explain differences in the results?
>In particular, why coefficients the same (the first result table
>columns),
>but standard errors is not?
>With best regards,
>Sergey S.
>
>	[[alternative HTML version deleted]]
>
>______________________________________________
>R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>https://stat.ethz.ch/mailman/listinfo/r-help
>PLEASE do read the posting guide
>http://www.R-project.org/posting-guide.html
>and provide commented, minimal, self-contained, reproducible code.

	[[alternative HTML version deleted]]



More information about the R-help mailing list