Estimating Quasi-Poisson Regression with GLIMMIX in SAS

R-bloggers 2015-10-15

(This article was first published on Yet Another Blog in Statistical Computing » S+/R, and kindly contributed to R-bloggers)

When modeling the frequency measure in the operational risk with regressions, most modelers often prefer Poisson or Negative Binomial regressions as best practices in the industry. However, as an alternative approach, Quasi-Poisson regression provides a more flexible model estimation routine with at least two benefits. First of all, Quasi-Poisson regression is able to address both over-dispersion and under-dispersion by assuming that the variance is a function of the mean such that VAR(Y|X) = Theta * MEAN(Y|X), where Theta > 1 for the over-dispersion and Theta < 1 for the under-dispersion. Secondly, estimated coefficients with Quasi-Poisson regression are identical to the ones with Standard Poisson regression, which is considered the prevailing practice in the industry.

While Quasi-Poisson regression can be easily estimated with glm() in R language, its estimation in SAS is not very straight-forward. Luckily, with GLIMMIX procedure, we can estimate Quasi-Poisson regression by directly specifying the functional relationship between the variance and the mean and making no distributional assumption in the MODEL statement, as demonstrated below.

proc glimmix data = credit_count;  model MAJORDRG = AGE ACADMOS MINORDRG OWNRENT / link = log solution;  _variance_ = _mu_;  random _residual_;run;  /*              Model Information Data Set                     WORK.CREDIT_COUNTResponse Variable            MAJORDRG        Response Distribution        Unknown         Link Function                Log             Variance Function            _mu_             Variance Matrix              Diagonal        Estimation Technique         Quasi-LikelihoodDegrees of Freedom Method    Residual                       Fit Statistics -2 Log Quasi-Likelihood           19125.57Quasi-AIC  (smaller is better)    19135.57Quasi-AICC (smaller is better)    19135.58Quasi-BIC  (smaller is better)    19173.10Quasi-CAIC (smaller is better)    19178.10Quasi-HQIC (smaller is better)    19148.09Pearson Chi-Square                51932.87Pearson Chi-Square / DF               3.86                        Parameter Estimates                         StandardEffect       Estimate       Error       DF    t Value    Pr > |t| Intercept     -1.3793     0.08613    13439     -16.01      <.0001AGE           0.01039    0.002682    13439       3.88      0.0001ACADMOS      0.001532    0.000385    13439       3.98      <.0001MINORDRG       0.4611     0.01348    13439      34.22      <.0001OWNRENT       -0.1994     0.05568    13439      -3.58      0.0003Residual       3.8643           .        .        .         .   */

For the comparison purpose, we also estimated a Quasi-Poisson regression in R, showing completely identical statistical results.

summary(glm(MAJORDRG ~ AGE + ACADMOS + MINORDRG + OWNRENT, data = credit_count, family = quasipoisson(link = "log")))  #               Estimate Std. Error t value Pr(>|t|)   # (Intercept) -1.3793249  0.0861324 -16.014  < 2e-16 ***# AGE          0.0103949  0.0026823   3.875 0.000107 ***# ACADMOS      0.0015322  0.0003847   3.983 6.84e-05 ***# MINORDRG     0.4611297  0.0134770  34.216  < 2e-16 ***# OWNRENT     -0.1993933  0.0556757  -3.581 0.000343 ***# ---# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1# # (Dispersion parameter for quasipoisson family taken to be 3.864409)# #     Null deviance: 24954  on 13443  degrees of freedom# Residual deviance: 22048  on 13439  degrees of freedom# AIC: NA

To leave a comment for the author, please follow the link and comment on their blog: Yet Another Blog in Statistical Computing » S+/R.

R-bloggers.com offers daily e-mail updates about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...