Formulas and counts now work!

Makefile

@@ -1,5 +1,5 @@
-update:
-	rsync -avz ../barry/include/barry include/
-
 build:
-	pip3 install .
+	pip3 install .
+
+update:
+	rsync -avz ../barry/include/barry include/

README.md

@@ -36,24 +36,57 @@ produces binary “wheels” for all platforms is illustrated in the
 import pydefm as m
 import numpy as np
 
-y = np.array([0, 10, 3])
-x = np.array([1, 2.0, 3.4])
-id = np.array([11, 2, 3])
+y = np.column_stack(
+  (np.array([0, 0, 1, 1, 1, 1]),
+  np.array([0, 1, 1, 0, 0, 1]))
+)
 
-obj = m.new_defm(id, y, x)
+x = np.array([1, 2.0, 3.4, 3, 1, 2])
+id = np.array([1, 1, 1, 2, 2, 2])
+
+obj = m.new_defm(id, y, x, column_major = False)
 
 # Printing the object on screen shows it is a pointer
 obj
 ```
 
-    <pydefm._core.DEFM at 0x7feeb40c6b70>
+    <pydefm._core.DEFM at 0x7fb88ff22070>
+
+Adding terms via formula
 
 ``` python
-# We can export member functions from the C++ class. Here is
-# an example of a function defined in the C++ class:
+m.term_formula(obj, "{y0}")
+m.term_formula(obj, "{y1}")
+m.term_formula(obj, "{0y0, y1}")
+obj.init()
 obj.print()
 ```
 
+``` python
+counts = m.get_stats(obj)
+counts
+```
+
+    array([[0., 0., 0.],
+           [0., 1., 1.],
+           [1., 1., 0.],
+           [1., 0., 0.],
+           [1., 0., 0.],
+           [1., 1., 0.]])
+
+Notice the first two columns coincide with the `y` array:
+
+``` python
+y
+```
+
+    array([[0, 0],
+           [0, 1],
+           [1, 1],
+           [1, 0],
+           [1, 0],
+           [1, 1]])
+
 Currently, the C++ code uses `printf` to print to the screen, which is a
 different buffer from Python (the problem will be solved using <a
 href="https://pybind11.readthedocs.io/en/stable/advanced/pycpp/utilities.html?highlight=print#using-python-s-print-function-in-c"
@@ -64,13 +97,17 @@ target="_blank">this</a>).
 m.print_y(obj)
 ```
 
-    Num. of Arrays       : 0
-    Support size         : 0
-    Support size range   : [2147483647, 0]
+    Num. of Arrays       : 6
+    Support size         : 6
+    Support size range   : [4, 4]
     Transform. Fun.      : no
-    Model terms (0)    :
-    Model Y variables (3):
+    Model terms (3)    :
+     - Motif {y0⁺}
+     - Motif {y1⁺}
+     - Motif {y0⁻, y1⁺}
+    Model rules (1)     :
+     - Markov model of order 0
+    Model Y variables (2):
       0) y0
       1) y1
-      2) y2
-    0 10 3 
+    0 0 

README.qmd

@@ -65,7 +65,7 @@ y = np.column_stack(
 x = np.array([1, 2.0, 3.4, 3, 1, 2])
 id = np.array([1, 1, 1, 2, 2, 2])
 
-obj = m.new_defm(id, y, x)
+obj = m.new_defm(id, y, x, column_major = False)
 
 # Printing the object on screen shows it is a pointer
 obj
@@ -75,15 +75,22 @@ obj
 Adding terms via formula
 
 ```{python}
+m.term_formula(obj, "{y0}")
+m.term_formula(obj, "{y1}")
 m.term_formula(obj, "{0y0, y1}")
 obj.init()
+obj.print() # Right now, printing shows at the end of the screen
 ```
 
 ```{python}
-# We can export member functions from the C++ class. Here is
-# an example of a function defined in the C++ class:
-obj.print()
-m.get_stats(obj)
+counts = m.get_stats(obj)
+counts
+```
+
+Notice the first two columns coincide with the `y` array:
+
+```{python}
+y
 ```
 
 Currently, the C++ code uses `printf` to print to the screen, which is a different buffer from Python (the problem will be solved using [this](https://pybind11.readthedocs.io/en/stable/advanced/pycpp/utilities.html?highlight=print#using-python-s-print-function-in-c){target="_blank"}).

include/barry/models/defm/defm-bones.hpp

@@ -87,6 +87,8 @@ class DEFM : public DEFMModel {
 
     std::vector< bool > is_motif();
 
+    bool get_column_major() const noexcept;
+
 };
 
 #endif

include/barry/models/defm/defm-meat.hpp

@@ -111,7 +111,7 @@ inline DEFM::DEFM(
     size_t m_order,
     bool copy_data,
     bool column_major
-) {
+) : column_major(column_major) {
 
     // Pointers
     if (copy_data)
@@ -445,6 +445,11 @@ inline std::vector< bool > DEFM::is_motif()
     return res;
 }
 
+inline bool DEFM::get_column_major() const noexcept
+{
+    return column_major;
+}
+
 #undef DEFM_RANGES
 #undef DEFM_LOOP_ARRAYS
 

src/model.cpp

@@ -35,25 +35,44 @@ namespace py = pybind11;
 py::array_t<double> get_stats(std::shared_ptr< defm::DEFM > m)
 {
 
-    auto model = m->get_model();
-
     // Getting sizes
-    size_t nrows = ptr->get_n_rows();
-    size_t ncols = model.nterms();
-    size_t m_ord = ptr->get_m_order();
+    size_t nrows = m->get_n_rows();
+    size_t ncols = m->nterms();
+    size_t m_ord = m->get_m_order();
 
-    const int * ID = ptr->get_ID();
+    const int * ID = m->get_ID();
 
     py::array_t< double > res({nrows, ncols});
-    auto target = model.get_stats_target();
 
+    auto res_buff = res.request();
+    double * res_ptr = static_cast< double * >(res_buff.ptr);
+
+    auto target = m->get_stats_target();
+
+    std::function<size_t(size_t,size_t,size_t,size_t)> element_access;
+
+    if (m->get_column_major())
+    {
+
+        element_access = [](size_t i, size_t j, size_t nrow, size_t) -> size_t {
+            return i + j * nrow;
+        };
 
+    } else {
+
+        element_access = [](size_t i, size_t j, size_t, size_t ncol) -> size_t {
+            return j + i * ncol;
+        };
+
+    }
+
+
     size_t i_effective = 0u;
     size_t n_obs_i = 0u;
     for (size_t i = 0u; i < nrows; ++i)
     {
 
-      // Do we need to reset the counter?
+        // Do we need to reset the counter?
         if ((i > 0) && (*(ID + i - 1u) != *(ID + i)))
             n_obs_i = 0u;
 
@@ -62,12 +81,17 @@ py::array_t<double> get_stats(std::shared_ptr< defm::DEFM > m)
         {
             // std::fill(res.row(i).begin(), res.row(i).end(), -99.0);
             for (size_t j = 0u; j < ncols; ++j)
-                res(i, j) = -99.0;
+            {
+                *(res_ptr + element_access(i, j, nrows, ncols)) = -99.0;
+                // res(i, j) = -99.0;
+            }
             continue;
         }
 
         for (size_t j = 0u; j < ncols; ++j)
-            res(i, j) = (*target)[i_effective][j];
+        {
+            *(res_ptr + element_access(i, j, nrows, ncols)) = (*target)[i_effective][j];
+        }
 
         i_effective++;