diff --git a/statsmodels/genmod/tests/results/results_glm.py b/statsmodels/genmod/tests/results/results_glm.py
index 7a5410edfe6..0b2a44035a2 100644
--- a/statsmodels/genmod/tests/results/results_glm.py
+++ b/statsmodels/genmod/tests/results/results_glm.py
@@ -85,6 +85,9 @@ def __init__(self):
         # TODO: taken from Stata; not available in sm yet
         self.chi2 = 1981.711859508729
 
+        self.prsquared = 0.90
+        self.prsquared_cox_snell = 1.00
+
         # self.pearson_chi2 = 836424.1293162981   # from Stata (?)
         self.fittedvalues = np.array([
             60055.659970240202, 61216.013942398131,
@@ -626,6 +629,8 @@ def __init__(self):
         self.llf_Stata = -2998.612927807218
         self.scale = 1.
         self.pearson_chi2 = 4051.921614
+        self.prsquared = 0.8346
+        self.prsquared_cox_snell = 1.0000
         self.resids = glm_test_resids.star98_resids
         self.fittedvalues = np.array([
             0.5833118,   0.75144661,  0.50058272, 0.68534524,  0.32251021,
@@ -795,6 +800,8 @@ def __init__(self):
         # self.llf = -82.47352  # Very close to ours as is
         self.scale = 0.003584283
         self.pearson_chi2 = .0860228056
+        self.prsquared = 0.429
+        self.prsquared_cox_snell = 0.97971
         self.resids = glm_test_resids.scotvote_resids
         self.fittedvalues = np.array([
             57.80431482,  53.2733447, 50.56347993, 58.33003783,
diff --git a/statsmodels/genmod/tests/test_glm.py b/statsmodels/genmod/tests/test_glm.py
index 857360115f1..5aa8da473be 100644
--- a/statsmodels/genmod/tests/test_glm.py
+++ b/statsmodels/genmod/tests/test_glm.py
@@ -191,6 +191,16 @@ def test_pearson_chi2(self):
             assert_allclose(self.res1.pearson_chi2, self.res2.pearson_chi2,
                             atol=1e-6, rtol=1e-6)
 
+    decimal_prsquared = DECIMAL_2
+    def test_prsquared(self):
+        if hasattr(self.res2, 'prsquared'):
+            assert_almost_equal(self.res1.prsquared, self.res2.prsquared, self.decimal_prsquared)
+
+    decimal_prsquared_cox_snell = DECIMAL_2
+    def test_prsquared_cox_snell(self):
+        if hasattr(self.res2, 'prsquared_cox_snell'):
+            assert_almost_equal(float(self.res1.prsquared_cox_snell), self.res2.prsquared_cox_snell, self.decimal_prsquared_cox_snell)
+
     @pytest.mark.smoke
     def test_summary(self):
         self.res1.summary()