fix: bug when fit_intercept=False, support.size=0

abess-team · Sep 5, 2023 · 79367b5 · 79367b5
2 parents 96ed869 + 5130593
commit 79367b5
Show file tree

Hide file tree

Showing 3 changed files with 23 additions and 6 deletions.
diff --git a/.github/workflows/r_website.yml b/.github/workflows/r_website.yml
@@ -35,9 +35,9 @@ jobs:
     steps:
       - uses: actions/checkout@v2
 
-      - uses: r-lib/actions/setup-r@v1
+      - uses: r-lib/actions/setup-r@v2
 
-      - uses: r-lib/actions/setup-pandoc@v1
+      - uses: r-lib/actions/setup-pandoc@v2
 
       - name: Query dependencies
         run: |

diff --git a/R-package/vignettes/v05-coxreg.Rmd b/R-package/vignettes/v05-coxreg.Rmd
@@ -160,6 +160,7 @@ fitted.results <- predict(abess_fit, newx = test, type = 'response')
 We now calculate the C-index, i.e., the probability that, for a pair of randomly chosen comparable samples, the sample with the higher risk prediction will experience an event before the other sample or belong to a higher binary class. On this dataset, the C-index is 0.64.
 ```{r}
 library(Hmisc)
+library(survival)
 Cindex <- max(1-rcorr.cens(fitted.results, Surv(test[, 1], test[,2]))[1],rcorr.cens(fitted.results, Surv(test[, 1], test[,2]))[1])
 Cindex
 ```

diff --git a/src/AlgorithmGLM.h b/src/AlgorithmGLM.h
@@ -336,7 +336,7 @@ class abessLogistic : public _abessGLM<Eigen::VectorXd, Eigen::VectorXd, double,
     }
 
     bool null_model(Eigen::VectorXd &y, Eigen::VectorXd &weights, double &coef0) {
-        coef0 = -log(1 / y.mean() - 1);
+        coef0 = -log(weights.sum() / y.dot(weights) - 1);
         return true;
     }
 };
@@ -406,7 +406,7 @@ class abessLm : public _abessGLM<Eigen::VectorXd, Eigen::VectorXd, double, T4> {
                            double loss0, Eigen::VectorXi &A, Eigen::VectorXi &g_index, Eigen::VectorXi &g_size) {
         // int n = x.rows();
         // int p = x.cols();
-        //if (x.cols() == 0) return true;
+        if (x.cols() == 0) return null_model(y, weights, coef0);
 
         // to ensure
         T4 X_full;
@@ -552,6 +552,11 @@ class abessLm : public _abessGLM<Eigen::VectorXd, Eigen::VectorXd, double, T4> {
             return enp;
         }
     }
+
+    bool null_model(Eigen::VectorXd &y, Eigen::VectorXd &weights, double &coef0) {
+        coef0 = y.dot(weights) / weights.sum();
+        return true;
+    }
 };
 
 template <class T4>
@@ -982,7 +987,7 @@ class abessMLm : public _abessGLM<Eigen::MatrixXd, Eigen::MatrixXd, Eigen::Vecto
         int n = x.rows();
         int p = x.cols();
         int M = y.cols();
-        //if (p == 0) return true;
+        if (p == 0) return null_model(y, weights, coef0);
 
         // to ensure
         T4 X;
@@ -1143,6 +1148,11 @@ class abessMLm : public _abessGLM<Eigen::MatrixXd, Eigen::MatrixXd, Eigen::Vecto
             return enp;
         }
     }
+
+    bool null_model(Eigen::MatrixXd &y, Eigen::VectorXd &weights, Eigen::VectorXd &coef0) {
+        coef0 = weights.transpose() * y / weights.sum();
+        return true;
+    }
 };
 
 template <class T4>
@@ -1171,7 +1181,7 @@ class abessMultinomial : public _abessGLM<Eigen::MatrixXd, Eigen::MatrixXd, Eige
         int n = x.rows();
         int p = x.cols();
         int M = y.cols();
-        //if (p == 0) return true;
+        if (p == 0) return null_model(y, weights, coef0);
 
         T4 X;
         add_constant_column(X, x, this->fit_intercept);
@@ -1511,6 +1521,12 @@ class abessMultinomial : public _abessGLM<Eigen::MatrixXd, Eigen::MatrixXd, Eige
             return enp;
         }
     }
+
+    bool null_model(Eigen::MatrixXd &y, Eigen::VectorXd &weights, Eigen::VectorXd &coef0) {
+        coef0 = (weights.transpose() * y / weights.sum()).array().log();
+        coef0 = coef0.array() - coef0(coef0.size() - 1);
+        return true;
+    }
 };
 
 template <class T4>