Adding class. accuracy as goodness of fit

python/logistic/logit_reg_boots.ipynb

@@ -587,7 +587,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 38,
+   "execution_count": null,
    "id": "2a860703",
    "metadata": {},
    "outputs": [
@@ -608,7 +608,8 @@
     "print(f\"{bpars_nonpar.shape = }\")\n",
     "\n",
     "# select a windows of parameters with the largest peak in the distribution\n",
-    "mask_sel = data_utils.within_bounds(bpars_nonpar, np.array([[-25,0], [0, 10]]))\n",
+    "window =  np.array([[-50,0], [0, 20]]\n",
+    "mask_sel = data_utils.within_bounds(bpars_nonpar,window))\n",
     "bpars_nonpar_sel = bpars_nonpar[mask_sel]\n",
     "\n",
     "print(\"size_selected:\", len(bpars_nonpar_sel))"
@@ -785,7 +786,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 43,
+   "execution_count": null,
    "id": "c53978ca",
    "metadata": {},
    "outputs": [
@@ -799,7 +800,8 @@
    ],
    "source": [
     "# select a windows of parameters with the largest peak in the distribution\n",
-    "mask_sel = data_utils.within_bounds(bpars_nonpar_strat, np.array([[-25,0], [0, 10]]))\n",
+    "window = np.array([[-100,0], [0, 20]]\n",
+    "mask_sel = data_utils.within_bounds(bpars_nonpar_strat, window)\n",
     "bpars_nonpar_strat_sel = bpars_nonpar_strat[mask_sel]\n",
     "\n",
     "print(\"size_selected:\", len(bpars_nonpar_strat_sel))"
@@ -910,7 +912,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 47,
+   "execution_count": null,
    "id": "c157ea73",
    "metadata": {},
    "outputs": [
@@ -924,7 +926,8 @@
    ],
    "source": [
     "# select a windows of parameters with the largest peak in the distribution\n",
-    "mask_sel = data_utils.within_bounds(bpars_param, np.array([[-25,0], [0, 10]]))\n",
+    "window = np.array([[-100,0], [0, 20]]\n",
+    "mask_sel = data_utils.within_bounds(bpars_param, window))\n",
     "bpars_param_sel = bpars_param[mask_sel]\n",
     "\n",
     "print(\"size_selected:\", len(bpars_param_sel))"

python/logistic/logit_utils.py

@@ -75,16 +75,21 @@ def logit_poly_fit(lm, x, y, degree = 1):
         x: array of n floats
         y: array of n int in {0,1}
         b: array of r = degree+1 floats, model parameters
-    
+        thresh: float, default 0.5, threshold value for classification
+
     Return:
-        {"n": n, "k":k, "dof":n-k, "LLF": log_likelihood, "AIC": AIC, "BIC": BIC}
+        {"n": n, "k":k, "dof":n-k, 
+        "LLF": log_likelihood, 
+        "AIC": AIC, 
+        "BIC": BIC, 
+        "A": classification accuracy (threshold values = 0.5 prob)}
     
     Ref:
         https://en.wikipedia.org/wiki/Logistic_regression
         https://en.wikipedia.org/wiki/Akaike_information_criterion
         https://www.medicine.mcgill.ca/epidemiology/joseph/courses/epib-621/logfit.pdf
 """
-def logit_poly_goodness_of_fit(x, y, b):
+def logit_poly_goodness_of_fit(x, y, b, thresh=0.5):
     
     # model probabilities
     p = logit_poly_model(x, b)
@@ -104,8 +109,12 @@ def logit_poly_goodness_of_fit(x, y, b):
     chi2 = np.sum(r**2)
     p_val = scipy.stats.chi2.sf(chi2, dof)
 
+    # using model as classifier
+    matches = y == np.heaviside(p - thresh, 1)
+
     return {"LLF": llf, "AIC": AIC, "BIC": BIC, 
-            "chi2": chi2, "p-value(chi2)": p_val,
+            "A" : np.count_nonzero (matches)/n,
+            "chi2": chi2, "p-value(chi2)": p_val,  # not very useful
             "n": n, "k": k, "dof": dof}
 
 """
@@ -249,7 +258,10 @@ def logit_poly_model_quantiles_delta(x, probs, mean_b, cov):
 
 """
     Generate m parameters via non-parametric bootstrapping with a minimal constraint 
-    that both groups should be present in the sampled data.
+    that both groups should be present in the sampled data:
+
+    boostrapped sample = (xb, yb)  by sampling with replacement pairs (x_i, y_i) 
+                                   with condition that yb can not be just 0 or just 1 
 
      Input:
         lm: linear_model.LogisticRegression
@@ -289,7 +301,12 @@ def get_nonparam_boots_pars(lm, x, y, m, degree = 1, seed = 1977):
 
 
 """
-    Generate m parameters via non-parametric stratified bootstrapping.
+    Generate m parameters via non-parametric stratified bootstrapping:
+
+        boostrapped sample = (xb, yb)  
+            
+            xb = (sampled with replacement from x0, sampled with replacement from x1)
+            yb = (0...0, 1...1)
 
     Input: 
         lm: linear_model.LogisticRegression
@@ -327,7 +344,11 @@ def get_nonparam_stratified_boots_pars(lm, x, y, m, degree = 1, seed = 1977):
     return np.array(lst)
 
 """
-    Generate m parameters via parametric bootstrapping.
+    Generate m parameters via parametric bootstrapping:
+
+        boostrapped sample = (x, yb) yb ~ B(model(x, fitted pars))  
+    
+    where B is Bernoulli distribution
 
     Input: 
         lm: linear_model.LogisticRegression
@@ -343,6 +364,7 @@ def get_nonparam_stratified_boots_pars(lm, x, y, m, degree = 1, seed = 1977):
     Ref:
 
       * https://www.scirp.org/journal/paperinformation?paperid=70962
+      * https://en.wikipedia.org/wiki/Bernoulli_distribution
 """
 def get_parametric_boots_pars(lm, x, y, m, degree = 1, seed = 1977):