diff --git a/src/Makefile.am b/src/Makefile.am
index 50a07b651bac9999a01a826977b28a3d7429e592..7ddaa8a0a483662cdb02e33d768be3a19e7135eb 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -35,6 +35,7 @@ __top_builddir__sisso_cpp_SOURCES = \
     inputs/InputParser.cpp \
     descriptor_identifier/Model/Model.cpp \
     descriptor_identifier/SISSORegressor.cpp \
+    python/bindings.cpp \
     main.cpp
 
 
diff --git a/src/descriptor_identifier/Model/Model.cpp b/src/descriptor_identifier/Model/Model.cpp
index 712c35302bd3417b29283e655b90c9a18a050618..268b07051f1e9522d5ebd108716159870af80b22 100644
--- a/src/descriptor_identifier/Model/Model.cpp
+++ b/src/descriptor_identifier/Model/Model.cpp
@@ -75,7 +75,7 @@ std::string Model::toString() const
         unit_rep << " + a" << std::to_string(ff) << " * " << _feats[ff]->expr();
     return unit_rep.str();
 }
-
+//
 std::ostream& operator<< (std::ostream& outStream, const Model& model)
 {
     outStream << model.toString();
@@ -166,3 +166,31 @@ void Model::test_to_file(std::string filename, std::vector<int> test_inds)
     }
     out_file_stream.close();
 }
+
+void Model::register_python()
+{
+    using namespace boost::python;
+    class_<Model>("Model", init<std::vector<double>, std::vector<double>, std::vector<model_node_ptr>, std::vector<int>, std::vector<int>>())
+        .def("predict", &Model::predict)
+        .def("fit", &Model::predict_train)
+        // .def(str(self))
+        .def_readonly("_n_samp_train", &Model::_n_samp_train)
+        .def_readonly("_n_samp_test", &Model::_n_samp_test)
+        .def_readonly("_n_dim", &Model::_n_dim)
+        .def_readonly("_feats", &Model::_feats)
+        .def_readonly("_coefs", &Model::_coefs)
+        .def_readonly("_prop_train", &Model::_prop_train)
+        .def_readonly("_prop_test", &Model::_prop_test)
+        .def_readonly("_train_error", &Model::_train_error)
+        .def_readonly("_test_error", &Model::_test_error)
+        .def_readonly("_D_train", &Model::_D_train)
+        .def_readonly("_D_test", &Model::_D_test)
+        .def_readonly("_prop_train_est", &Model::_prop_train_est)
+        .def_readonly("_prop_test_est", &Model::_prop_test_est)
+        .def_readonly("_task_sizes_train", &Model::_task_sizes_train)
+        .def_readonly("_task_sizes_test", &Model::_task_sizes_test)
+        .add_property("rmse", &Model::rmse)
+        .add_property("test_rmse", &Model::test_rmse)
+        .add_property("max_ae", &Model::max_ae)
+        .add_property("test_max_ae", &Model::test_max_ae);
+}
\ No newline at end of file
diff --git a/src/descriptor_identifier/Model/Model.hpp b/src/descriptor_identifier/Model/Model.hpp
index 0c5ecb83e50c2e2782354c5d8626a2346467c53f..7f13c1c3a50146ac23a6e05a3539067fb936611f 100644
--- a/src/descriptor_identifier/Model/Model.hpp
+++ b/src/descriptor_identifier/Model/Model.hpp
@@ -2,6 +2,7 @@
 #define MODEL
 
 #include <boost/filesystem.hpp>
+#include <boost/python.hpp>
 
 #include<iomanip>
 #include<fstream>
@@ -55,12 +56,12 @@ public:
     /**
      * @brief Accessor function to _prop_est
      */
-    inline std::vector<double>& predict(){return _prop_test_est;}
+    inline std::vector<double> predict(){return _prop_test_est;}
 
     /**
      * @brief Accessor function to _prop_est
      */
-    inline std::vector<double>& predict_train(){return _prop_train_est;}
+    inline std::vector<double> predict_train(){return _prop_train_est;}
 
     /**
      * @brief Copy the error into a new array
@@ -98,6 +99,8 @@ public:
      * @brief Print model to a file
      */
     void train_to_file(std::string filename);
+
+    static void register_python();
 };
 
 /**
diff --git a/src/descriptor_identifier/SISSORegressor.cpp b/src/descriptor_identifier/SISSORegressor.cpp
index 21482fede7a6cd67ccb343243e6e6b364119bbb0..03d4c3b6df4809cca6fe55dbc41804905e584e71 100644
--- a/src/descriptor_identifier/SISSORegressor.cpp
+++ b/src/descriptor_identifier/SISSORegressor.cpp
@@ -3,12 +3,10 @@
 SISSORegressor::SISSORegressor(std::shared_ptr<FeatureSpace> feat_space, std::vector<double> prop, std::vector<double> prop_test, std::vector<int> task_sizes_train, std::vector<int> task_sizes_test, int n_dim, int n_residual):
     _prop(prop),
     _prop_test(prop_test),
-    _a(new double[(n_dim + 1) * prop.size()]),
-    _b(new double[prop.size()]),
-    _ones(new double[prop.size()]),
-    _error(new double[prop.size()]),
-    _work(nullptr),
-    _s(new double[n_dim + 1]),
+    _a((n_dim + 1) * prop.size()),
+    _b(prop.size()),
+    _error(prop.size(), 0.0),
+    _s(n_dim + 1),
     _task_sizes_train(task_sizes_train),
     _task_sizes_test(task_sizes_test),
     _feat_space(feat_space),
@@ -21,25 +19,20 @@ SISSORegressor::SISSORegressor(std::shared_ptr<FeatureSpace> feat_space, std::ve
 {
 
     // Initialize a, b, ones, s, and _error arrays
-    std::vector<double> zeros((_n_dim + 1) * _n_samp, 0.0);
-    std::copy_n(zeros.begin(), (_n_dim + 1) * _n_samp, _a.get());
-    std::copy_n(zeros.begin(), _n_samp, _b.get());
-    std::copy_n(zeros.begin(), _n_samp, _error.get());
-    std::copy_n(zeros.begin(), _n_dim + 1, _s.get());
-
-    std::vector<double> ones(_n_samp, 1.0);
-    std::copy_n(ones.begin(), ones.size(), _ones.get());
+    std::fill_n(_a.data(), (_n_dim + 1) * _n_samp, 0.0);
+    std::fill_n(_b.data(), _n_samp, 0.0);
+    std::fill_n(_s.data(), _n_dim + 1, 0.0);
 
     // // Get optimal work size
     _lwork = get_opt_lwork(_n_dim + 1);
-    _work = std::unique_ptr<double[]>(new double[_lwork]);
+    _work = std::vector<double>(_lwork, 0.0);
 }
 
 void  SISSORegressor::set_a(std::vector<int>& inds, int start, int n_samp)
 {
     for(int ii = 0; ii < inds.size(); ++ii)
-        std::copy_n(node_value_arrs::get_d_matrix_ptr(inds[ii]) + start, n_samp, _a.get() + ii * n_samp);
-    std::copy_n(_ones.get(), n_samp, _a.get() + inds.size() * n_samp);
+        std::copy_n(node_value_arrs::get_d_matrix_ptr(inds[ii]) + start, n_samp, _a.data() + ii * n_samp);
+    std::fill_n(_a.data() + inds.size() * n_samp, n_samp, 1.0);
 }
 
 void SISSORegressor::least_squares(std::vector<int>& inds, double* coeffs, int start, int n_samp)
@@ -48,13 +41,13 @@ void SISSORegressor::least_squares(std::vector<int>& inds, double* coeffs, int s
     int n_dim = inds.size() + 1;
 
     set_a(inds, start, n_samp);
-    std::copy_n(_prop.data() + start, n_samp, _b.get());
+    std::copy_n(_prop.data() + start, n_samp, _b.data());
 
-    dgelss_(n_samp, n_dim, 1, _a.get(), n_samp, _b.get(), n_samp, _s.get(), 1e-13, &_rank, _work.get(), _lwork, &info);
+    dgelss_(n_samp, n_dim, 1, _a.data(), n_samp, _b.data(), n_samp, _s.data(), 1e-13, &_rank, _work.data(), _lwork, &info);
 
     if(info == 0)
     {
-        std::copy_n(_b.get(), n_dim, coeffs);
+        std::copy_n(_b.data(), n_dim, coeffs);
     }
     else
     {
@@ -70,7 +63,7 @@ int SISSORegressor::get_opt_lwork(int n_dim)
     std::vector<double> work(1, 0.0);
     int info = 0, rank = 0;
 
-    dgelss_(_n_samp, n_dim, 1, _a.get(), _n_samp, _b.get(), _n_samp, _s.get(), 0.0, &rank, work.data(), -1, &info);
+    dgelss_(_n_samp, n_dim, 1, _a.data(), _n_samp, _b.data(), _n_samp, _s.data(), 0.0, &rank, work.data(), -1, &info);
 
     if(info == 0)
         return static_cast<int>(std::round(work[0]));
@@ -80,11 +73,11 @@ int SISSORegressor::get_opt_lwork(int n_dim)
 
 void SISSORegressor::set_error(std::vector<int>& inds, double* coeffs, int start, int n_samp)
 {
-    std::copy_n(_prop.begin() + start, n_samp, _error.get() + start);
+    std::fill_n(_error.data() + start, n_samp, -1.0 * coeffs[inds.size()]);
+    daxpy_(n_samp, 1.0,  _prop.data() + start, 1, _error.data() + start, 1);
 
     for(int ii = 0; ii < inds.size(); ++ii)
-        daxpy_(n_samp, -1.0 * coeffs[ii], node_value_arrs::get_d_matrix_ptr(inds[ii]) + start, 1, _error.get() + start, 1);
-    daxpy_(n_samp, -1.0 * coeffs[inds.size()], _ones.get(), 1, _error.get() + start, 1);
+        daxpy_(n_samp, -1.0 * coeffs[ii], node_value_arrs::get_d_matrix_ptr(inds[ii]) + start, 1, _error.data() + start, 1);
 }
 
 void SISSORegressor::fit()
@@ -161,8 +154,8 @@ void SISSORegressor::l0_norm(std::vector<double>& prop, int n_dim)
         {
             least_squares(inds, coefs.data(), start, sz);
             set_error(inds, coefs.data(), start, sz);
-            error += ddot_(sz, _error.get() + start, 1, _error.get() + start, 1) / sz;
-            // error += std::pow(util_funcs::norm(_error.get() + start, sz), 2.0) / sz;
+            error += ddot_(sz, _error.data() + start, 1, _error.data() + start, 1) / sz;
+            // error += std::pow(util_funcs::norm(_error.data() + start, sz), 2.0) / sz;
             start += sz;
         }
         error = std::sqrt(error / _task_sizes_train.size());
@@ -200,3 +193,21 @@ void SISSORegressor::l0_norm(std::vector<double>& prop, int n_dim)
 
     _models.push_back(models);
 }
+
+void SISSORegressor::register_python()
+{
+    using namespace boost::python;
+    class_<SISSORegressor>("SISSORegressor", init<std::shared_ptr<FeatureSpace>, std::vector<double>, std::vector<double>, std::vector<int>, std::vector<int>, int, int>())
+        .def("fit", &SISSORegressor::fit)
+        .add_property("prop", &SISSORegressor::prop)
+        .add_property("prop_test", &SISSORegressor::prop_test)
+        .add_property("models", &SISSORegressor::models)
+        .add_property("n_samp", &SISSORegressor::n_samp)
+        .add_property("n_dim", &SISSORegressor::n_dim)
+        .add_property("n_residual", &SISSORegressor::n_residual)
+        .add_property("feat_space", &SISSORegressor::feat_space)
+        .add_property("error", &SISSORegressor::error)
+        .def_readonly("_task_sizes_train", &SISSORegressor::_task_sizes_train)
+        .def_readonly("_task_sizes_test", &SISSORegressor::_task_sizes_test)
+        .def_readonly("_mpi_comm", &SISSORegressor::_mpi_comm);
+}
\ No newline at end of file
diff --git a/src/descriptor_identifier/SISSORegressor.hpp b/src/descriptor_identifier/SISSORegressor.hpp
index 1c87f3e467aa83b9fc460a951ed5e5f60091f865..f47bdf1f860dca5ed358200c8b0fbd39713d6384 100644
--- a/src/descriptor_identifier/SISSORegressor.hpp
+++ b/src/descriptor_identifier/SISSORegressor.hpp
@@ -16,13 +16,12 @@ protected:
 
     std::vector<double> _prop; //!< Property array
     std::vector<double> _prop_test; //!< Property array
+    std::vector<double> _error; //!< Array to calculate the residuals for the models
 
-    std::unique_ptr<double[]> _a; //!< A matrix for least squares
-    std::unique_ptr<double[]> _b; //!< Solution array for least squares
-    std::unique_ptr<double[]> _ones; //!< Array of ones to copy over for least squares comparison
-    std::unique_ptr<double[]> _error; //!< Array to calculate the residuals for the models
-    std::unique_ptr<double[]> _work; //!< The work array for least squares problems
-    std::unique_ptr<double[]> _s; //!< The S array for least squares problems
+    std::vector<double> _a; //!< A matrix for least squares
+    std::vector<double> _b; //!< Solution array for least squares
+    std::vector<double> _work; //!< The work array for least squares problems
+    std::vector<double> _s; //!< The S array for least squares problems
 
     std::vector<int> _task_sizes_train;
     std::vector<int> _task_sizes_test;
@@ -98,12 +97,17 @@ public:
     /**
      * @brief Acessor function for prop
      */
-    inline std::vector<double>& prop(){return _prop;}
+    inline std::vector<double> prop(){return _prop;}
+
+    /**
+     * @brief Acessor function for prop
+     */
+    inline std::vector<double> prop_test(){return _prop_test;}
 
     /**
      * @brief Acessor function for models
      */
-    inline std::vector<std::vector<Model>>& models(){return _models;}
+    inline std::vector<std::vector<Model>> models(){return _models;}
 
     /**
      * @brief Acessor function for n_samp
@@ -120,21 +124,12 @@ public:
      */
     inline int n_residual(){return _n_residual;}
 
-    /**
-     * @brief Acessor function for lwork
-     */
-    inline int lwork(){return _lwork;}
-
-    /**
-     * @brief Acessor function for rank
-     */
-    inline int rank(){return _rank;}
-
     /**
      * @brief Acessor function for {
      */
-    inline double* error(){return _error.get();}
+    inline std::vector<double> error(){return _error;}
 
+    static void register_python();
 };
 
 #endif
diff --git a/src/feature_creation/feature_space/FeatureSpace.cpp b/src/feature_creation/feature_space/FeatureSpace.cpp
index 8b673b96c6362911d3889e035aafbb7e7ed2082c..b4282d5e8a4ead534596de06824d128f25135266 100644
--- a/src/feature_creation/feature_space/FeatureSpace.cpp
+++ b/src/feature_creation/feature_space/FeatureSpace.cpp
@@ -383,14 +383,14 @@ void FeatureSpace::project_generated(double* prop, int size, std::vector<node_pt
             {
                 if(scores_sel.size() == _n_sis_select)
                 {
-                    phi_sel[worst_score_ind]->selected() = false;
-                    phi_sel[worst_score_ind]->set_dmat_ind(-1);
+                    phi_sel[worst_score_ind]->set_selected(false);
+                    phi_sel[worst_score_ind]->set_d_mat_ind(-1);
 
                     phi_sel[worst_score_ind] = generated_phi[inds[ii]];
                     scores_sel[worst_score_ind] = cur_score;
 
-                    phi_sel[worst_score_ind]->selected() = true;
-                    phi_sel[worst_score_ind]->set_dmat_ind(node_value_arrs::N_SELECTED - _n_sis_select + worst_score_ind);
+                    phi_sel[worst_score_ind]->set_selected(true);
+                    phi_sel[worst_score_ind]->set_d_mat_ind(node_value_arrs::N_SELECTED - _n_sis_select + worst_score_ind);
                     phi_sel[worst_score_ind]->set_value();
                 }
                 else
@@ -399,8 +399,8 @@ void FeatureSpace::project_generated(double* prop, int size, std::vector<node_pt
                     phi_sel.push_back(generated_phi[inds[ii]]);
                     scores_sel.push_back(cur_score);
 
-                    phi_sel.back()->selected() = true;
-                    phi_sel.back()->set_dmat_ind(node_value_arrs::N_SELECTED - _n_sis_select + scores_sel.size());
+                    phi_sel.back()->set_selected(true);
+                    phi_sel.back()->set_d_mat_ind(node_value_arrs::N_SELECTED - _n_sis_select + scores_sel.size());
                     phi_sel.back()->set_value();
                 }
                 worst_score_ind = std::max_element(scores_sel.begin(), scores_sel.end()) - scores_sel.begin();
@@ -486,8 +486,8 @@ void FeatureSpace::sis(std::vector<double>& prop)
             scores_sel[cur_feat_local] = _scores[inds[ii]];
             // phi_sel.push_back(std::make_shared<FeatureNode>(cur_feat + cur_feat_local, _phi[inds[ii]]->expr(), _phi[inds[ii]]->value(), _phi[inds[ii]]->test_value(), _phi[inds[ii]]->unit(), true));
             phi_sel.push_back(_phi[inds[ii]]);
-            phi_sel.back()->selected() = true;
-            phi_sel.back()->set_dmat_ind(cur_feat + cur_feat_local);
+            phi_sel.back()->set_selected(true);
+            phi_sel.back()->set_d_mat_ind(cur_feat + cur_feat_local);
             phi_sel.back()->set_value();
             ++cur_feat_local;
         }
@@ -553,18 +553,18 @@ void FeatureSpace::sis(std::vector<double>& prop)
             if(std::none_of(inds.begin(), inds.end(), [&compare_ind](int i1){return i1 == compare_ind;}))
             {
                 scores_sel.erase(scores_sel.begin() + ii);
-                phi_sel[ii]->selected() = false;
-                phi_sel[ii]->set_dmat_ind(-1);
+                phi_sel[ii]->set_selected(false);
+                phi_sel[ii]->set_d_mat_ind(-1);
                 phi_sel.erase(phi_sel.begin() + ii);
             }
             else
             {
-                phi_sel[ii]->selected() = false;
-                phi_sel[ii]->set_dmat_ind(-1);
+                phi_sel[ii]->set_selected(false);
+                phi_sel[ii]->set_d_mat_ind(-1);
             }
         }
         for(auto& feat : phi_sel)
-            feat->selected() = false;
+            feat->set_selected(false);
 
 
         scores_sel.resize(_n_sis_select, 0.0);
@@ -592,8 +592,8 @@ void FeatureSpace::sis(std::vector<double>& prop)
                 {
                     out_file_stream << std::setw(feat_num_width) <<std::left << cur_feat << std::setw(24) << std::setprecision(18) << std::left << sent_scores[inds[ii]] << sent_phi[inds[ii]]->expr() << std::endl;
                     _phi_selected.push_back(sent_phi[inds[ii]]);
-                    _phi_selected.back()->selected() = true;
-                    _phi_selected.back()->set_dmat_ind(cur_feat);
+                    _phi_selected.back()->set_selected(true);
+                    _phi_selected.back()->set_d_mat_ind(cur_feat);
                     ++cur_feat;
                 }
             }
@@ -617,8 +617,8 @@ void FeatureSpace::sis(std::vector<double>& prop)
 
                         _phi_selected.push_back(sent_phi[inds[ii]]);
 
-                        _phi_selected.back()->selected() = true;
-                        _phi_selected.back()->set_dmat_ind(cur_feat);
+                        _phi_selected.back()->set_selected(true);
+                        _phi_selected.back()->set_d_mat_ind(cur_feat);
                         _phi_selected.back()->set_value();
 
                         scores_sel[cur_feat_local] = sent_scores[inds[ii]];
@@ -660,7 +660,7 @@ void FeatureSpace::sis(std::vector<double>& prop)
         {
             out_file_stream << std::setw(feat_num_width) <<std::left << cur_feat << std::setw(24) << std::setprecision(18) << std::left << scores_sel[ind] << phi_sel[ind]->expr() << std::endl;
             _phi_selected.push_back(phi_sel[ind]);
-            _phi_selected.back()->set_dmat_ind(cur_feat);
+            _phi_selected.back()->set_d_mat_ind(cur_feat);
             _phi_selected.back()->set_value();
             ++cur_feat;
             ++cur_feat_local;
@@ -676,3 +676,44 @@ void FeatureSpace::sis(std::vector<double>& prop)
         out_file_stream << std::endl;
     }
 }
+
+void FeatureSpace::register_python()
+{
+    using namespace boost::python;
+    class_<FeatureSpace>("FeatureSpace", init<std::shared_ptr<MPI_Interface>, std::vector<node_ptr>, std::vector<std::string>, std::vector<double>, std::vector<int>, optional<int, int, int, int, int, double, double>>())
+        .def("sis", &FeatureSpace::sis)
+        .def("feat_in_phi", &FeatureSpace::feat_in_phi)
+        .add_property("phi_selected", &FeatureSpace::phi_selected)
+        .add_property("phi", &FeatureSpace::phi)
+        .add_property("phi0", &FeatureSpace::phi0)
+        .add_property("scores", &FeatureSpace::scores)
+        .add_property("task_sizes", &FeatureSpace::task_sizes)
+        .def_readonly("_allowed_ops", &FeatureSpace::_allowed_ops)
+        .def_readonly("_un_operators", &FeatureSpace::_un_operators)
+        .def_readonly("_com_bin_operators", &FeatureSpace::_com_bin_operators)
+        .def_readonly("_bin_operators", &FeatureSpace::_bin_operators)
+        .def_readonly("_start_gen", &FeatureSpace::_start_gen)
+        .def_readonly("_feature_space_file", &FeatureSpace::_feature_space_file)
+        .def_readonly("_l_bound", &FeatureSpace::_l_bound)
+        .def_readonly("_u_bound", &FeatureSpace::_u_bound)
+        .def_readonly("_max_phi", &FeatureSpace::_max_phi)
+        .def_readonly("_n_sis_select", &FeatureSpace::_n_sis_select)
+        .def_readonly("_n_samp", &FeatureSpace::_n_samp)
+        .def_readonly("_n_feat", &FeatureSpace::_n_feat)
+        .def_readonly("_n_rung_store", &FeatureSpace::_n_rung_store)
+        .def_readonly("_n_rung_generate", &FeatureSpace::_n_rung_generate)
+        .def_readonly("_mpi_comm", &FeatureSpace::_mpi_comm);
+}
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/src/feature_creation/feature_space/FeatureSpace.hpp b/src/feature_creation/feature_space/FeatureSpace.hpp
index d59fab2532588312bd58fdf58155e94226c77f24..a741a565ccbdf256dc395fd6ad75e005bbf6d27a 100644
--- a/src/feature_creation/feature_space/FeatureSpace.hpp
+++ b/src/feature_creation/feature_space/FeatureSpace.hpp
@@ -9,6 +9,7 @@
 
 #include <boost/serialization/shared_ptr.hpp>
 #include <boost/filesystem.hpp>
+#include <boost/python.hpp>
 
 #include <iostream>
 #include <iomanip>
@@ -100,7 +101,7 @@ public:
     /**
      * @brief Accessor function for _scores
      */
-    inline std::vector<double>& scores(){return _scores;};
+    inline std::vector<double> scores(){return _scores;};
 
     /**
      * @brief Accessor function for _mpi_comm
@@ -132,6 +133,7 @@ public:
      */
     inline bool feat_in_phi(int ind){return (ind >= _phi[0]->feat_ind()) && (ind <= _phi.back()->feat_ind());}
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/FeatureNode.cpp b/src/feature_creation/node/FeatureNode.cpp
index 1428c4ece12f21e406b9b84ef71ee5292c68225a..457ee7f54267c76bc896ad029d6122c1a3a0cea2 100644
--- a/src/feature_creation/node/FeatureNode.cpp
+++ b/src/feature_creation/node/FeatureNode.cpp
@@ -37,4 +37,24 @@ void FeatureNode::update_div_mult_leaves(std::map<std::string, double>& div_mult
     expected_abs_tot += std::abs(fact);
 }
 
+void FeatureNode::register_python()
+{
+    std::string (FeatureNode::*expr_1)() = &FeatureNode::expr;
+    std::string (FeatureNode::*expr_const)() const = &FeatureNode::expr;
+
+    using namespace boost::python;
+    class_<FeatureNode, bases<Node>>("FeatureNode", init<int, std::string, std::vector<double>, std::vector<double>, Unit>())
+        .def("is_nan", &FeatureNode::is_nan)
+        .def("is_const", &FeatureNode::is_const)
+        .def("set_value", &FeatureNode::set_value)
+        .def("set_test_value", &FeatureNode::set_test_value)
+        .add_property("value", &FeatureNode::value)
+        .add_property("test_value", &FeatureNode::test_value)
+        .add_property("expr", expr_1)
+        .add_property("expr", expr_const)
+        .add_property("unit", &FeatureNode::unit)
+        .add_property("type", &FeatureNode::type)
+        .add_property("rung", &FeatureNode::rung)
+    ;
+}
 // BOOST_CLASS_EXPORT(FeatureNode)
diff --git a/src/feature_creation/node/FeatureNode.hpp b/src/feature_creation/node/FeatureNode.hpp
index 4f7bef882d28daa0e63cf3fed89c18d60e8842b9..0eaaa89da54a1ba4df1005e9a3cf53713da9754c 100644
--- a/src/feature_creation/node/FeatureNode.hpp
+++ b/src/feature_creation/node/FeatureNode.hpp
@@ -146,6 +146,8 @@ public:
      * @param add_sub_leaves [description]
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
+
+    static void register_python();
 };
 
 #endif
diff --git a/src/feature_creation/node/ModelNode.cpp b/src/feature_creation/node/ModelNode.cpp
index dfe69e7dc3bcc8c7457af23b441b8dd0609db748..cc48ac80f0c92e016e98884b14d2da8377f04847 100644
--- a/src/feature_creation/node/ModelNode.cpp
+++ b/src/feature_creation/node/ModelNode.cpp
@@ -6,7 +6,7 @@ ModelNode::ModelNode()
 ModelNode::ModelNode(int feat_ind, int rung, std::string expr, std::vector<double> value, std::vector<double> test_value, Unit unit) :
     FeatureNode(feat_ind, expr, value, test_value, unit),
     _rung(rung)
-{}
+ {}
 
 ModelNode::~ModelNode()
 {}
@@ -31,4 +31,18 @@ void ModelNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_l
     expected_abs_tot += std::abs(fact);
 }
 
-// BOOST_CLASS_EXPORT(ModelNode)
+void ModelNode::register_python()
+{
+    std::string (ModelNode::*expr_1)() = &ModelNode::expr;
+    std::string (ModelNode::*expr_const)() const = &ModelNode::expr;
+
+    using namespace boost::python;
+    class_<ModelNode, bases<FeatureNode>>("ModelNode", init<int, int, std::string, std::vector<double>, std::vector<double>, Unit>())
+        .def("is_nan", &ModelNode::is_nan)
+        .def("is_const", &ModelNode::is_const)
+        .def("set_value", &ModelNode::set_value)
+        .def("set_test_value", &ModelNode::set_test_value)
+        .add_property("type", &ModelNode::type)
+        .add_property("rung", &ModelNode::rung)
+    ;
+}
diff --git a/src/feature_creation/node/ModelNode.hpp b/src/feature_creation/node/ModelNode.hpp
index d80a2f8ee09e9d784edb054c30afa26266f8681d..c945fa4efade4e3c6e7ef1a9ba1ccd3c3ce74a53 100644
--- a/src/feature_creation/node/ModelNode.hpp
+++ b/src/feature_creation/node/ModelNode.hpp
@@ -108,6 +108,8 @@ public:
      * @param add_sub_leaves [description]
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
+
+    static void register_python();
 };
 
 #endif
diff --git a/src/feature_creation/node/Node.cpp b/src/feature_creation/node/Node.cpp
index 22432f60db9b30290c3e5c6317b2f5b5a0cca688..c11fabee24032a5d570934399307225f3e660f16 100644
--- a/src/feature_creation/node/Node.cpp
+++ b/src/feature_creation/node/Node.cpp
@@ -15,5 +15,106 @@ Node::Node(int feat_ind, int n_samp, int n_test_samp) :
 Node::~Node()
 {}
 
+using namespace boost::python;
+struct NodeWrap : Node, wrapper<Node>
+{
+    std::string expr()
+    {
+        return this->get_override("expr")();
+    }
+
+    Unit unit()
+    {
+        return this->get_override("unit")();
+    }
+
+    std::vector<double> value()
+    {
+        return this->get_override("value")();
+    }
+
+    std::vector<double> test_value()
+    {
+        return this->get_override("test_value")();
+    }
+
+    void set_value(int offset = -1)
+    {
+        this->get_override("set_value")();
+    }
+
+    double* value_ptr(int offset = -1)
+    {
+        return this->get_override("value_ptr")();
+    }
+
+    void set_test_value(int offset = -1)
+    {
+        this->get_override("set_test_value")();
+    }
+
+    double* test_value_ptr(int offset = -1)
+    {
+        return this->get_override("test_value_ptr")();
+    }
+
+    bool is_nan()
+    {
+        return this->get_override("is_nan")();
+    }
+
+    bool is_const()
+    {
+        return this->get_override("is_const")();
+    }
+
+    NODE_TYPE type()
+    {
+        return this->get_override("type")();
+    }
+
+    int rung(int cur_rung = 0)
+    {
+        return this->get_override("rung")();
+    }
+
+    void update_add_sub_leaves(std::map<std::string, int>& add_sub_leaves, int pl_mn, int& expected_abs_tot)
+    {
+        this->get_override("update_add_sub_leaves");
+    }
+
+    void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot)
+    {
+        this->get_override("update_div_mult_leaves");
+    }
+
+};
+
+void Node::register_python()
+{
+    void (Node::*reindex_1)(int) = &Node::reindex;
+    void (Node::*reindex_2)(int, int) = &Node::reindex;
+    class_<NodeWrap, boost::noncopyable>("Node", no_init)
+        .def("reindex", reindex_1)
+        .def("reindex", reindex_2)
+        .add_property("n_samp", &Node::n_samp)
+        .add_property("n_test_samp", &Node::n_test_samp)
+        .add_property("feat_ind", &Node::feat_ind)
+        .add_property("arr_ind", &Node::arr_ind)
+        .add_property("selected", &Node::selected, &Node::set_selected)
+        .add_property("d_mat_ind", &Node::d_mat_ind, &Node::set_d_mat_ind)
+        .def("expr", pure_virtual(&Node::expr))
+        .def("unit", pure_virtual(&Node::unit))
+        .def("value", pure_virtual(&Node::value))
+        .def("test_value", pure_virtual(&Node::test_value))
+        .def("set_value", pure_virtual(&Node::set_value))
+        .def("set_test_value", pure_virtual(&Node::set_test_value))
+        .def("is_nan", pure_virtual(&Node::is_nan))
+        .def("is_const", pure_virtual(&Node::is_const))
+        .def("type", pure_virtual(&Node::type))
+        .def("rung", pure_virtual(&Node::rung))
+        ;
+}
+
 BOOST_SERIALIZATION_ASSUME_ABSTRACT(Node)
 
diff --git a/src/feature_creation/node/Node.hpp b/src/feature_creation/node/Node.hpp
index 5e28fa08703402e476dde69884751419b3c39e10..131dd1591024e124fab34780a9aa9f135019ac29 100644
--- a/src/feature_creation/node/Node.hpp
+++ b/src/feature_creation/node/Node.hpp
@@ -101,16 +101,18 @@ public:
     /**
      * @brief Accessor function to get the feature ind
      */
-    inline int& feat_ind(){return _feat_ind;}
+    inline int feat_ind(){return _feat_ind;}
 
     /**
      * @brief Accessor function to get the feature array index
      */
-    inline int& arr_ind(){return _arr_ind;}
+    inline int arr_ind(){return _arr_ind;}
 
-    inline bool& selected(){return _selected;}
+    inline bool selected(){return _selected;}
 
-    inline void set_dmat_ind(int ind){_d_mat_ind = ind;}
+    inline void set_selected(bool sel){_selected = sel;}
+
+    inline void set_d_mat_ind(int ind){_d_mat_ind = ind;}
 
     inline int d_mat_ind(){return _d_mat_ind;}
     /**
@@ -190,6 +192,8 @@ public:
      * @param expected_abs_tot The expected absolute sum of all values in div_mult_leaves
      */
     virtual void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot) = 0;
+
+    static void register_python();
 };
 
 typedef std::shared_ptr<Node> node_ptr;
diff --git a/src/feature_creation/node/operator_nodes/OperatorNode.hpp b/src/feature_creation/node/operator_nodes/OperatorNode.hpp
index 6e358ceae62c49bd0d6dfcb69cdb2be602aaa867..e8f4a578e88cce496305f06adc0e27026ce31699 100644
--- a/src/feature_creation/node/operator_nodes/OperatorNode.hpp
+++ b/src/feature_creation/node/operator_nodes/OperatorNode.hpp
@@ -12,6 +12,7 @@
 #include <boost/serialization/array.hpp>
 #include <boost/serialization/vector.hpp>
 
+using namespace boost::python;
 
 /**
  * @brief Base class to describe operator nodes
@@ -166,6 +167,66 @@ public:
      */
     virtual void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot) = 0;
 
+    struct OperatorNodeWrap : OperatorNode<N>, wrapper<OperatorNode<N>>
+    {
+
+        void set_value(int offset = -1)
+        {
+            this->get_override("set_value")();
+        }
+
+        void set_test_value(int offset = -1)
+        {
+            this->get_override("set_test_value")();
+        }
+
+        NODE_TYPE type()
+        {
+            return this->get_override("type")();
+        }
+
+        int rung(int cur_rung = 0)
+        {
+            return this->get_override("rung")();
+        }
+
+        std::string expr()
+        {
+            return this->get_override("expr")();
+        }
+
+        Unit unit()
+        {
+            return this->get_override("unit")();
+        }
+
+        void update_add_sub_leaves(std::map<std::string, int>& add_sub_leaves, int pl_mn, int& expected_abs_tot)
+        {
+            this->get_override("update_add_sub_leaves")();
+        }
+        void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot)
+        {
+            this->get_override("update_div_mult_leaves")();
+        }
+    };
+
+    static void register_python()
+    {
+        // class_<OperatorNode, bases<Node>>("OperatorNode", init<std::array<node_ptr, N>, int>())
+        // ;
+        class_<OperatorNodeWrap, bases<Node>, boost::noncopyable>("OperatorNode")
+            .def("is_nan", &OperatorNode::is_nan)
+            .def("is_const", &OperatorNode::is_const)
+            .add_property("value", &OperatorNode::value)
+            .add_property("test_value", &OperatorNode::test_value)
+            .def("set_value", pure_virtual(&OperatorNode::set_value))
+            .def("set_test_value", pure_virtual(&OperatorNode::set_test_value))
+            .def("type", pure_virtual(&OperatorNode::type))
+            .def("rung", pure_virtual(&OperatorNode::rung))
+            .def("expr", pure_virtual(&OperatorNode::expr))
+            .def("unit", pure_virtual(&OperatorNode::unit))
+        ;
+    }
 };
 
 #endif
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.cpp
index a201993c8f1add6b8c62510a773ee753e067a368..59a62cfd0a4ed7af31f48d224825772cc1773282 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.cpp
@@ -32,7 +32,7 @@ AbsDiffNode::AbsDiffNode(std::array<node_ptr, 2> feats, int feat_ind, double l_b
         throw InvalidFeatureException();
 
     set_test_value();
- }
+}
 
 AbsDiffNode::AbsDiffNode(node_ptr feat_1, node_ptr feat_2, int feat_ind, double l_bound, double u_bound):
     OperatorNode({feat_1, feat_2}, feat_ind)
@@ -84,3 +84,17 @@ void AbsDiffNode::update_div_mult_leaves(std::map<std::string, double>& div_mult
 
     expected_abs_tot += std::abs(fact);
 }
+
+void AbsDiffNode::register_python()
+{
+    using namespace boost::python;
+    class_<AbsDiffNode, bases<OperatorNode<2>>>("AbsDiffNode", init<node_ptr, node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 2>, int, double, double>())
+        .def("set_value", &AbsDiffNode::set_value)
+        .def("set_test_value", &AbsDiffNode::set_test_value)
+        .add_property("expr", &AbsDiffNode::expr)
+        .add_property("unit", &AbsDiffNode::unit)
+        .add_property("type", &AbsDiffNode::type)
+        .add_property("rung", &AbsDiffNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.hpp
index 994e80298d9d422f37c8ff3d53aee2b6d1087f80..4eab084620c715a8437b15611a504f5bd1a88b4f 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_difference.hpp
@@ -67,6 +67,8 @@ public:
      * @param expected_abs_tot The expected absolute sum of all values in div_mult_leaves
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
+
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.cpp
index 90ee52bf8f6c2bec4591159db37807a1fa4ebd1d..fe476da4d42ae4c5bb1e63376b93cb534e71c081 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.cpp
@@ -46,3 +46,17 @@ void AbsNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
 
     expected_abs_tot += std::abs(fact);
 }
+
+void AbsNode::register_python()
+{
+    using namespace boost::python;
+    class_<AbsNode, bases<OperatorNode<1>>>("AbsNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &AbsNode::set_value)
+        .def("set_test_value", &AbsNode::set_test_value)
+        .add_property("expr", &AbsNode::expr)
+        .add_property("unit", &AbsNode::unit)
+        .add_property("type", &AbsNode::type)
+        .add_property("rung", &AbsNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.hpp
index 44d7fc2a189ddcb4126acb7118e0e9ff8fcf71f2..80479324fbe354c385747cd9fe6010a24a10724e 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/absolute_value.hpp
@@ -66,6 +66,7 @@ public:
      * @param expected_abs_tot The expected absolute sum of all values in div_mult_leaves
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.cpp
index 4786123b978114a1143d1edc17a6683759be4743..d621608432e4cc6c185425cd182fb94286be6578 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.cpp
@@ -29,7 +29,7 @@ AddNode::AddNode(std::array<node_ptr, 2> feats, int feat_ind, double l_bound, do
         throw InvalidFeatureException();
 
     set_test_value();
- }
+}
 
 AddNode::AddNode(node_ptr feat_1, node_ptr feat_2, int feat_ind, double l_bound, double u_bound):
     OperatorNode({feat_1, feat_2}, feat_ind)
@@ -57,7 +57,7 @@ AddNode::AddNode(node_ptr feat_1, node_ptr feat_2, int feat_ind, double l_bound,
         throw InvalidFeatureException();
 
     set_test_value();
- }
+}
 
 void AddNode::update_add_sub_leaves(std::map<std::string, int>& add_sub_leaves, int pl_mn, int& expected_abs_tot)
 {
@@ -75,3 +75,17 @@ void AddNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
 
     expected_abs_tot += std::abs(fact);
 }
+
+void AddNode::register_python()
+{
+    using namespace boost::python;
+    class_<AddNode, bases<OperatorNode<2>>>("AddNode", init<node_ptr, node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 2>, int, double, double>())
+        .def("set_value", &AddNode::set_value)
+        .def("set_test_value", &AddNode::set_test_value)
+        .add_property("expr", &AddNode::expr)
+        .add_property("unit", &AddNode::unit)
+        .add_property("type", &AddNode::type)
+        .add_property("rung", &AddNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.hpp
index f00fe1c4f127a37fd3d92c69f81e31b113dacd90..04296e3c7595c9a8690d3f7a28cd1994de56d59d 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add.hpp
@@ -67,6 +67,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.cpp
index 9bf1fab672a151297a12d09282f9d7bb0bd01d7e..2e14154afa7c0283323578dc4d5a15b3d21d57d7 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.cpp
@@ -17,7 +17,7 @@ CosNode::CosNode(std::array<node_ptr, 1> feats, int feat_ind, double l_bound, do
         throw InvalidFeatureException();
 
        set_test_value();
- }
+}
 
 CosNode::CosNode(node_ptr feat, int feat_ind, double l_bound, double u_bound):
     OperatorNode({feat}, feat_ind)
@@ -33,7 +33,7 @@ CosNode::CosNode(node_ptr feat, int feat_ind, double l_bound, double u_bound):
         throw InvalidFeatureException();
 
        set_test_value();
- }
+}
 
 void CosNode::update_add_sub_leaves(std::map<std::string, int>& add_sub_leaves, int pl_mn, int& expected_abs_tot)
 {
@@ -57,3 +57,16 @@ void CosNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     expected_abs_tot += std::abs(fact);
 }
 
+void CosNode::register_python()
+{
+    using namespace boost::python;
+    class_<CosNode, bases<OperatorNode<1>>>("CosNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &CosNode::set_value)
+        .def("set_test_value", &CosNode::set_test_value)
+        .add_property("expr", &CosNode::expr)
+        .add_property("unit", &CosNode::unit)
+        .add_property("type", &CosNode::type)
+        .add_property("rung", &CosNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.hpp
index b40670a9f2b0487fb5a7ae832f77245e8dbba397..c701887866e5b6a94adc3d6c4c0837f84c0733bd 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos.hpp
@@ -67,6 +67,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.cpp
index 6d53cea31bed572543f1a7d9951616cb9b9eae84..055b3361a51d82951900b1a504832aad7b2fbc47 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.cpp
@@ -44,3 +44,17 @@ void CbNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_leav
 {
     _feats[0]->update_div_mult_leaves(div_mult_leaves, fact * 3.0, expected_abs_tot);
 }
+
+void CbNode::register_python()
+{
+    using namespace boost::python;
+    class_<CbNode, bases<OperatorNode<1>>>("CbNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &CbNode::set_value)
+        .def("set_test_value", &CbNode::set_test_value)
+        .add_property("expr", &CbNode::expr)
+        .add_property("unit", &CbNode::unit)
+        .add_property("type", &CbNode::type)
+        .add_property("rung", &CbNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.hpp
index 5dd0071b6ec93e3135ea89a636023c252d501543..d1e25bb1bbbcc0ab959256e944fc9c16e73b96c6 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube.hpp
@@ -66,6 +66,7 @@ public:
      * @param expected_abs_tot The expected absolute sum of all values in div_mult_leaves
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.cpp
index 47f4daff556b126c7fd163a0ed2b42ca5a1cdcdd..6cf69aeeadd744cb4ff665d41bca0c15367c9121 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.cpp
@@ -44,3 +44,17 @@ void CbrtNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_le
 {
     _feats[0]->update_div_mult_leaves(div_mult_leaves, fact / 3.0, expected_abs_tot);
 }
+
+void CbrtNode::register_python()
+{
+    using namespace boost::python;
+    class_<CbrtNode, bases<OperatorNode<1>>>("CbrtNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &CbrtNode::set_value)
+        .def("set_test_value", &CbrtNode::set_test_value)
+        .add_property("expr", &CbrtNode::expr)
+        .add_property("unit", &CbrtNode::unit)
+        .add_property("type", &CbrtNode::type)
+        .add_property("rung", &CbrtNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.hpp
index 5c2c636e942f0c25dc2b5abccf59a9679b6a1f11..5d6ed1621ebcb24ffbe163859ae2a473232c2d3b 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cube_root.hpp
@@ -67,6 +67,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.cpp
index 098e623cf20a60d54bcbfd992c45fe31e6bc05e0..d07c20c2282483e3170ee11dd58b10f937d6e9c5 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.cpp
@@ -75,3 +75,17 @@ void DivNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     _feats[0]->update_div_mult_leaves(div_mult_leaves, fact, expected_abs_tot);
     _feats[1]->update_div_mult_leaves(div_mult_leaves, -1.0*fact, expected_abs_tot);
 }
+
+void DivNode::register_python()
+{
+    using namespace boost::python;
+    class_<DivNode, bases<OperatorNode<2>>>("DivNode", init<node_ptr, node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 2>, int, double, double>())
+        .def("set_value", &DivNode::set_value)
+        .def("set_test_value", &DivNode::set_test_value)
+        .add_property("expr", &DivNode::expr)
+        .add_property("unit", &DivNode::unit)
+        .add_property("type", &DivNode::type)
+        .add_property("rung", &DivNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.hpp
index 769a17393e40c1b6dcf841091aabff6138365321..9a6cabb95b6682d54b136da8c91526abcb6eadf0 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/divide.hpp
@@ -67,6 +67,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.cpp
index 13f4f5b7f0bbf598e752f2b7f774628233c1ebea..5d97e99b3abfa99b45611a524311a9be945325ad 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.cpp
@@ -17,7 +17,7 @@ ExpNode::ExpNode(std::array<node_ptr, 1> feats, int feat_ind, double l_bound, do
         throw InvalidFeatureException();
 
        set_test_value();
- }
+}
 
 ExpNode::ExpNode(node_ptr feat, int feat_ind, double l_bound, double u_bound):
     OperatorNode({feat}, feat_ind)
@@ -33,7 +33,7 @@ ExpNode::ExpNode(node_ptr feat, int feat_ind, double l_bound, double u_bound):
         throw InvalidFeatureException();
 
        set_test_value();
- }
+}
 
 void ExpNode::update_add_sub_leaves(std::map<std::string, int>& add_sub_leaves, int pl_mn, int& expected_abs_tot)
 {
@@ -56,3 +56,17 @@ void ExpNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
 
     expected_abs_tot += std::abs(fact);
 }
+
+void ExpNode::register_python()
+{
+    using namespace boost::python;
+    class_<ExpNode, bases<OperatorNode<1>>>("ExpNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &ExpNode::set_value)
+        .def("set_test_value", &ExpNode::set_test_value)
+        .add_property("expr", &ExpNode::expr)
+        .add_property("unit", &ExpNode::unit)
+        .add_property("type", &ExpNode::type)
+        .add_property("rung", &ExpNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.hpp
index 138dde1421dfd946af1a2c3b97b0b11c87d27518..49f4bdb9bff0c6bbcf786e51ad4d936bb08b6c29 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exponential.hpp
@@ -67,6 +67,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.cpp
index 12d2d7d3df7a2713545097cc5378bb08b439320a..5e35fdb6a3dce90c75cf3aaff09ea7f747788573 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.cpp
@@ -44,3 +44,17 @@ void InvNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
 {
     _feats[0]->update_div_mult_leaves(div_mult_leaves, fact * -1.0, expected_abs_tot);
 }
+
+void InvNode::register_python()
+{
+    using namespace boost::python;
+    class_<InvNode, bases<OperatorNode<1>>>("InvNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &InvNode::set_value)
+        .def("set_test_value", &InvNode::set_test_value)
+        .add_property("expr", &InvNode::expr)
+        .add_property("unit", &InvNode::unit)
+        .add_property("type", &InvNode::type)
+        .add_property("rung", &InvNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.hpp
index bc942c01b2789d4a8659181b2daf400fb666286c..9303c75bf015612afa7d6bed6f1e448d558f6c2e 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inverse.hpp
@@ -68,5 +68,6 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.cpp
index 112709b2fbd91ee01378c1d2bfc3daa6284b9dba..7b6eab0cb4289403c890cb13e712fb8ad3322eb4 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.cpp
@@ -56,3 +56,17 @@ void LogNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
 
     expected_abs_tot += std::abs(fact);
 }
+
+void LogNode::register_python()
+{
+    using namespace boost::python;
+    class_<LogNode, bases<OperatorNode<1>>>("LogNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &LogNode::set_value)
+        .def("set_test_value", &LogNode::set_test_value)
+        .add_property("expr", &LogNode::expr)
+        .add_property("unit", &LogNode::unit)
+        .add_property("type", &LogNode::type)
+        .add_property("rung", &LogNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.hpp
index 81b97ba7722051fd390d40ee5c26d94aee3590c4..af70a775bab2361174229dfd587fb90e8ad3e9d9 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log.hpp
@@ -66,6 +66,7 @@ public:
      * @param expected_abs_tot The expected absolute sum of all values in div_mult_leaves
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.cpp
index 16aed94d8da1652889283d04e320a3c1b917b62d..66a8bbd1f40f69c755d708cbc02a4b42b338d70c 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.cpp
@@ -75,3 +75,17 @@ void MultNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_le
     _feats[0]->update_div_mult_leaves(div_mult_leaves, fact, expected_abs_tot);
     _feats[1]->update_div_mult_leaves(div_mult_leaves, fact, expected_abs_tot);
 }
+
+void MultNode::register_python()
+{
+    using namespace boost::python;
+    class_<MultNode, bases<OperatorNode<2>>>("MultNode", init<node_ptr, node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 2>, int, double, double>())
+        .def("set_value", &MultNode::set_value)
+        .def("set_test_value", &MultNode::set_test_value)
+        .add_property("expr", &MultNode::expr)
+        .add_property("unit", &MultNode::unit)
+        .add_property("type", &MultNode::type)
+        .add_property("rung", &MultNode::rung)
+    ;
+}
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.hpp
index a937ca35233d5f65eae8b7547c5a0ad01f4cc0e4..a7be03d55543b635372fd93a5a7b58b363a14044 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/multiply.hpp
@@ -68,6 +68,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.cpp
index 063c5dbdd73d66f0b96e162201f39577cfc85d23..c92c426a6f18c3efd055bd27dff7d8dcf5535a51 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.cpp
@@ -56,3 +56,17 @@ void NegExpNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_
 
     expected_abs_tot += std::abs(fact);
 }
+
+void NegExpNode::register_python()
+{
+    using namespace boost::python;
+    class_<NegExpNode, bases<OperatorNode<1>>>("NegExpNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &NegExpNode::set_value)
+        .def("set_test_value", &NegExpNode::set_test_value)
+        .add_property("expr", &NegExpNode::expr)
+        .add_property("unit", &NegExpNode::unit)
+        .add_property("type", &NegExpNode::type)
+        .add_property("rung", &NegExpNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.hpp
index 81f622c02d9a9b9701783cef5bb01383fedf22f6..664f581b5c73ef1615affd6a4b6b7d15543da3cb 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/negative_exponential.hpp
@@ -68,6 +68,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.cpp
index 8bb02076c009281fba23b8e99991b5f4239d5641..02e1f41534142facb97b36dcfea63533377905b3 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.cpp
@@ -56,3 +56,17 @@ void SinNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
 
     expected_abs_tot += std::abs(fact);
 }
+
+void SinNode::register_python()
+{
+    using namespace boost::python;
+    class_<SinNode, bases<OperatorNode<1>>>("SinNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &SinNode::set_value)
+        .def("set_test_value", &SinNode::set_test_value)
+        .add_property("expr", &SinNode::expr)
+        .add_property("unit", &SinNode::unit)
+        .add_property("type", &SinNode::type)
+        .add_property("rung", &SinNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.hpp
index 1463ef34b815bd3f81a6c08b5ff21b4bc2122f40..275366af44b1b0c03cdf539857692f4d2e14ac95 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin.hpp
@@ -68,6 +68,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.cpp
index 2e67b4b7d9bd7cceea7a7d3dfbc6e2ef9d5f696a..c66bd122244599c088fd6da2d490b20455543ca6 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.cpp
@@ -44,3 +44,17 @@ void SixPowNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_
 {
     _feats[0]->update_div_mult_leaves(div_mult_leaves, fact * 6.0, expected_abs_tot);
 }
+
+void SixPowNode::register_python()
+{
+    using namespace boost::python;
+    class_<SixPowNode, bases<OperatorNode<1>>>("SixPowNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &SixPowNode::set_value)
+        .def("set_test_value", &SixPowNode::set_test_value)
+        .add_property("expr", &SixPowNode::expr)
+        .add_property("unit", &SixPowNode::unit)
+        .add_property("type", &SixPowNode::type)
+        .add_property("rung", &SixPowNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.hpp
index 770af9f940e6665cbdf485a442fb7a6fd535bea8..bad3ff5f7de7d76d3894f4122c2854da7e900c23 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sixth_power.hpp
@@ -68,6 +68,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.cpp
index f85e1f0bb547f5290304b0cdfb7db4b4741c5f3e..49a83ddc768cef7b0d1e531733237fc3e141fa37 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.cpp
@@ -44,3 +44,17 @@ void SqNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_leav
 {
     _feats[0]->update_div_mult_leaves(div_mult_leaves, fact * 2.0, expected_abs_tot);
 }
+
+void SqNode::register_python()
+{
+    using namespace boost::python;
+    class_<SqNode, bases<OperatorNode<1>>>("SqNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &SqNode::set_value)
+        .def("set_test_value", &SqNode::set_test_value)
+        .add_property("expr", &SqNode::expr)
+        .add_property("unit", &SqNode::unit)
+        .add_property("type", &SqNode::type)
+        .add_property("rung", &SqNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.hpp
index 7c46c86d5104d8c811b64f083ee73caa8c88977a..850acfd5e6edf9d740c5a6e522506ca5f167d86d 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square.hpp
@@ -66,6 +66,7 @@ public:
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
 
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.cpp
index ea54c2cd7578bce3439fa57109a3f239a32fe8fc..2a5a606c95226e00db1b916be2c4f4fab90bb6b4 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.cpp
@@ -44,3 +44,17 @@ void SqrtNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_le
 {
     _feats[0]->update_div_mult_leaves(div_mult_leaves, fact / 2.0, expected_abs_tot);
 }
+
+void SqrtNode::register_python()
+{
+    using namespace boost::python;
+    class_<SqrtNode, bases<OperatorNode<1>>>("SqrtNode", init<node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 1>, int, double, double>())
+        .def("set_value", &SqrtNode::set_value)
+        .def("set_test_value", &SqrtNode::set_test_value)
+        .add_property("expr", &SqrtNode::expr)
+        .add_property("unit", &SqrtNode::unit)
+        .add_property("type", &SqrtNode::type)
+        .add_property("rung", &SqrtNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.hpp
index 3f327423669998ce8b3ead9e4532b33f4924690f..e06b3bf01518b25426e7dc6a8f46a7ef2a2fdf60 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/square_root.hpp
@@ -67,6 +67,7 @@ public:
      * @param expected_abs_tot The expected absolute sum of all values in div_mult_leaves
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.cpp
index f3a654766d2d022cf2c711fdf39fb429b57e2ef9..2ed49e6c5c4e80fc9dd0a905a3ed489492ab8544 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.cpp
@@ -28,8 +28,8 @@ SubNode::SubNode(std::array<node_ptr, 2> feats, int feat_ind, double l_bound, do
     if(is_nan() || is_const() || (util_funcs::max_abs_val<double>(value_ptr(), _n_samp) > u_bound) || (util_funcs::max_abs_val<double>(value_ptr(), _n_samp) < l_bound))
         throw InvalidFeatureException();
 
-       set_test_value();
- }
+    set_test_value();
+}
 
 SubNode::SubNode(node_ptr feat_1, node_ptr feat_2, int feat_ind, double l_bound, double u_bound):
     OperatorNode({feat_1, feat_2}, feat_ind)
@@ -75,3 +75,17 @@ void SubNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
 
     expected_abs_tot += std::abs(fact);
 }
+
+void SubNode::register_python()
+{
+    using namespace boost::python;
+    class_<SubNode, bases<OperatorNode<2>>>("SubNode", init<node_ptr, node_ptr, int, double, double>())
+        .def(init<std::array<node_ptr, 2>, int, double, double>())
+        .def("set_value", &SubNode::set_value)
+        .def("set_test_value", &SubNode::set_test_value)
+        .add_property("expr", &SubNode::expr)
+        .add_property("unit", &SubNode::unit)
+        .add_property("type", &SubNode::type)
+        .add_property("rung", &SubNode::rung)
+    ;
+}
\ No newline at end of file
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.hpp
index 40b87a26f95f034473ab936701d74be108315938..49ee521d310e88a73880ed8d4d3658b2be7a379f 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/subtract.hpp
@@ -67,6 +67,7 @@ public:
      * @param expected_abs_tot The expected absolute sum of all values in div_mult_leaves
      */
     void update_div_mult_leaves(std::map<std::string, double>& div_mult_leaves, double fact, double& expected_abs_tot);
+    static void register_python();
 };
 
 #endif
\ No newline at end of file
diff --git a/src/feature_creation/units/Unit.cpp b/src/feature_creation/units/Unit.cpp
index c384ea8583f3591f784833413a286f3643f356e6..61637af2f279b08f24ae0989c1ef0b09afa56b46 100644
--- a/src/feature_creation/units/Unit.cpp
+++ b/src/feature_creation/units/Unit.cpp
@@ -157,3 +157,21 @@ std::ostream& operator<< (std::ostream& outStream, const Unit& unit)
     outStream << unit.toString();
     return outStream;
 }
+
+
+void Unit::register_python()
+{
+    using namespace boost::python;
+    class_<Unit>("Unit", init<>())
+        .def(init<std::map<std::string, double>>())
+        .def(init<std::string>())
+        .def(init<Unit&>())
+        // .def(str(self))
+        .def("inverse", &Unit::inverse)
+        .def(self * self)
+        .def(self / self)
+        .def(self == self)
+        .def(self != self)
+        // .def(pow(self, other<double>))
+        .add_property("dct", &Unit::dct);
+}
\ No newline at end of file
diff --git a/src/feature_creation/units/Unit.hpp b/src/feature_creation/units/Unit.hpp
index 9f360f9ffdf293014f92258b3a59d686f4f3a633..94f9074712e38303185a1c2aef6bf0edf7f823ed 100644
--- a/src/feature_creation/units/Unit.hpp
+++ b/src/feature_creation/units/Unit.hpp
@@ -12,6 +12,7 @@
 #include <boost/algorithm/string.hpp>
 #include <boost/serialization/map.hpp>
 #include <boost/serialization/string.hpp>
+#include <boost/python.hpp>
 
 using StringRange = boost::iterator_range<std::string::const_iterator>;
 
@@ -126,7 +127,10 @@ public:
         ar & _dct;
     }
 
+    static void register_python();
 };
+
+
 /**
  * @brief Operator to print the unit to a string stream
  *
diff --git a/src/mpi_interface/MPI_Interface.cpp b/src/mpi_interface/MPI_Interface.cpp
index 4c18b6f03fd54f0721235fca41634202d2cf932f..ddb8b988b4b008d4b1ea3089985db43044683f6d 100644
--- a/src/mpi_interface/MPI_Interface.cpp
+++ b/src/mpi_interface/MPI_Interface.cpp
@@ -14,3 +14,29 @@ std::array<int, 2> MPI_Interface::get_start_end_for_iterator(int length, int sta
 
     return start_end;
 }
+
+std::shared_ptr<MPI_Interface> mpi_setup::comm;
+
+void mpi_setup::init_mpi_env(int& argc, char **&argv)
+{
+    if(env == 0)
+        env = new boost::mpi::environment(argc, argv);
+}
+
+void mpi_setup::init_mpi_env()
+{
+    if(env == 0)
+    {
+        #ifdef BOOST_MPI_HAS_NOARG_INITIALIZATION
+            env = new boost::mpi::environment();
+        #else
+            throw std::runtime_error("MPI cannot be initialized without arguments");
+        #endif
+    }
+}
+
+void mpi_setup::finalize_mpi_env()
+{
+    delete env;
+    env = 0;
+}
\ No newline at end of file
diff --git a/src/mpi_interface/MPI_Interface.hpp b/src/mpi_interface/MPI_Interface.hpp
index 0024b28d3e109273c26696bcd6965af57c2c608e..329e92e7c019400420c8e0d6ba95568e6894d512 100644
--- a/src/mpi_interface/MPI_Interface.hpp
+++ b/src/mpi_interface/MPI_Interface.hpp
@@ -4,7 +4,6 @@
 #include <boost/mpi.hpp>
 #include <array>
 
-extern boost::mpi::environment env;
 namespace mpi = boost::mpi;
 
 // Augment the boost mpi communicator class with some other useful data
@@ -13,7 +12,7 @@ namespace mpi = boost::mpi;
  * @details MPI communicator used to transfer data throughout the cell.
  *
  */
-class MPI_Interface : public boost::mpi::communicator
+class MPI_Interface : public mpi::communicator
 {
 public:
     /**
@@ -38,4 +37,16 @@ public:
 
 };
 
+namespace mpi_setup
+{
+    static mpi::environment* env = 0;
+    extern std::shared_ptr<MPI_Interface> comm;
+
+    void init_mpi_env(int& argc, char **&argv);
+
+    void init_mpi_env();
+
+    void finalize_mpi_env();
+}
+
 #endif
\ No newline at end of file
diff --git a/src/python/_sisso.cpp b/src/python/_sisso.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..86653ab4783b8a77dca199f4953238e4db363b18
--- /dev/null
+++ b/src/python/_sisso.cpp
@@ -0,0 +1,17 @@
+#include <python/bindings.hpp>
+#include <utils/MPI_Interface.hpp>
+#include <Python.h>
+
+static void finalize();
+
+BOOST_PYTHON_MODULE(_sisso)
+{
+    mpi_setup::init_mpi_env();
+    sisso::register_python();
+    Py_AtExit(&finalize);
+}
+
+void finalize()
+{
+    mpi_setup::finalize_mpi_env();
+}
\ No newline at end of file
diff --git a/src/python/_sisso.hpp b/src/python/_sisso.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/src/python/bindings.cpp b/src/python/bindings.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..d79bdf65222fcdc21bbd30acf22e5f229d823930
--- /dev/null
+++ b/src/python/bindings.cpp
@@ -0,0 +1,31 @@
+#include <python/bindings.hpp>
+
+void sisso::register_python()
+{
+    Model::register_python();
+    SISSORegressor::register_python();
+    FeatureSpace::register_python();
+    Unit::register_python();
+    Node::register_python();
+    FeatureNode::register_python();
+    ModelNode::register_python();
+    OperatorNode<1>::register_python();
+    OperatorNode<2>::register_python();
+    AddNode::register_python();
+    SubNode::register_python();
+    DivNode::register_python();
+    MultNode::register_python();
+    AbsDiffNode::register_python();
+    AbsNode::register_python();
+    InvNode::register_python();
+    LogNode::register_python();
+    ExpNode::register_python();
+    NegExpNode::register_python();
+    SinNode::register_python();
+    CosNode::register_python();
+    CbNode::register_python();
+    CbrtNode::register_python();
+    SqNode::register_python();
+    SqrtNode::register_python();
+    SixPowNode::register_python();
+}
\ No newline at end of file
diff --git a/src/python/bindings.hpp b/src/python/bindings.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..ee80b8c6553c182cdab158cdcfc16d21339e9bd8
--- /dev/null
+++ b/src/python/bindings.hpp
@@ -0,0 +1,13 @@
+
+#ifndef PYTHON_BINDINGS
+#define PYTHON_BINDINGS
+
+#include <descriptor_identifier/SISSORegressor.hpp>
+#include <feature_creation/feature_space/FeatureSpace.hpp>
+
+namespace sisso
+{
+    void register_python();
+}
+
+#endif
\ No newline at end of file