diff --git a/src/feature_creation/feature_space/FeatureSpace.cpp b/src/feature_creation/feature_space/FeatureSpace.cpp
index 57a66471c1550bb16b89a6b4544bff6144bf74df..ff09826c38e5c3922dba1f77913b6e87a403eafe 100644
--- a/src/feature_creation/feature_space/FeatureSpace.cpp
+++ b/src/feature_creation/feature_space/FeatureSpace.cpp
@@ -234,36 +234,30 @@ void FeatureSpace::generate_feature_space(std::vector<double>& prop)
                 std::vector<int> del_inds;
 
                 _mpi_comm->barrier();
-                #pragma omp parallel
+                node_value_arrs::clear_temp_reg();
+                for(int sc = 0; sc < _scores.size() - 1; ++sc)
                 {
-                    std::vector<int> del_inds_private;
-                    #pragma omp for schedule(dynamic)
-                    for(int sc = 0; sc < _scores.size() - 1; ++sc)
+                    if(_scores[inds[sc]] > -1e-10)
                     {
-                        if(_scores[inds[sc]] > -1e-10)
+                        double base_val = std::abs(util_funcs::r(_phi[_start_gen.back() + inds[sc]]->value_ptr(), _phi[_start_gen.back() + inds[sc]]->value_ptr(), _n_samp));
+                        for(int sc2 = sc + 1; sc2 < _scores.size(); ++sc2)
                         {
-                            for(int sc2 = sc + 1; sc2 < _scores.size(); ++sc2)
+                            if(std::abs(base_val - std::abs(util_funcs::r(_phi[_start_gen.back() + inds[sc]]->value_ptr(), _phi[_start_gen.back() + inds[sc2]]->value_ptr(0, true), _n_samp))) < 1e-13)
                             {
-                                if(std::abs(util_funcs::r(_phi[_start_gen.back() + inds[sc]]->value_ptr(0), _phi[_start_gen.back() + inds[sc]]->value_ptr(1), _n_samp) - std::abs(util_funcs::r(_phi[_start_gen.back() + inds[sc]]->value_ptr(0), _phi[_start_gen.back() + inds[sc2]]->value_ptr(1), _n_samp))) < 1e-13)
-                                {
-                                    del_inds_private.push_back(-1 * (inds[sc] + _start_gen.back()));
-                                    break;
-                                }
+                                del_inds.push_back(-1 * (inds[sc] + _start_gen.back()));
+                                break;
                             }
                         }
-                        else if(_scores[inds[sc + 1]] - _scores[inds[sc]] < 1e-10)
-                        {
-                            if(std::abs(util_funcs::r(_phi[_start_gen.back() + inds[sc]]->value_ptr(0), _phi[_start_gen.back() + inds[sc]]->value_ptr(1), _n_samp) - std::abs(util_funcs::r(_phi[_start_gen.back() + inds[sc]]->value_ptr(0), _phi[_start_gen.back() + inds[sc + 1]]->value_ptr(1), _n_samp))) < 1e-13)
-                                del_inds_private.push_back(-1 * (inds[sc] + _start_gen.back()));
-                        }
                     }
-
-                    #pragma omp critical
-                    del_inds.insert(del_inds.end(), del_inds_private.begin(), del_inds_private.end());
+                    else if(_scores[inds[sc + 1]] - _scores[inds[sc]] < 1e-10)
+                    {
+                        double base_val = std::abs(util_funcs::r(_phi[_start_gen.back() + inds[sc]]->value_ptr(), _phi[_start_gen.back() + inds[sc]]->value_ptr(), _n_samp));
+                        if(std::abs(base_val - std::abs(util_funcs::r(_phi[_start_gen.back() + inds[sc]]->value_ptr(), _phi[_start_gen.back() + inds[sc + 1]]->value_ptr(0, true), _n_samp))) < 1e-13)
+                            del_inds.push_back(-1 * (inds[sc] + _start_gen.back()));
+                    }
                 }
 
                 inds = util_funcs::argsort(del_inds);
-
                 for(int ii = 0; ii < inds.size(); ++ii)
                     _phi.erase(_phi.begin() - del_inds[inds[ii]]);
 
@@ -271,6 +265,8 @@ void FeatureSpace::generate_feature_space(std::vector<double>& prop)
                 for(int ff = _start_gen.back(); ff < _phi.size(); ++ff)
                     _phi[ff]->reindex(ff);
             }
+            node_value_arrs::clear_temp_reg();
+
             if(nn <= _n_rung_store)
             {
                 node_value_arrs::resize_values_arr(nn, _phi.size(), (nn != _max_phi));
@@ -472,14 +468,14 @@ void FeatureSpace::project_generated(double* prop, int size, std::vector<node_pt
                     if(scores_sel_private.size() == _n_sis_select)
                     {
                         generated_phi[inds[ii]]->reindex(node_value_arrs::N_STORE_FEATURES + worst_score_ind + _n_sis_select * (omp_get_thread_num() + _mpi_comm->size()));
-                        // generated_phi[inds[ii]]->set_value();
+                        generated_phi[inds[ii]]->set_value();
                         phi_sel_private[worst_score_ind] = generated_phi[inds[ii]];
                         scores_sel_private[worst_score_ind] = cur_score;
                     }
                     else
                     {
                         generated_phi[inds[ii]]->reindex(node_value_arrs::N_STORE_FEATURES + scores_sel_private.size() + _n_sis_select * (omp_get_thread_num() +  _mpi_comm->size()));
-                        // generated_phi[inds[ii]]->set_value();
+                        generated_phi[inds[ii]]->set_value();
                         phi_sel_private.push_back(generated_phi[inds[ii]]);
                         scores_sel_private.push_back(cur_score);
                     }
@@ -494,7 +490,7 @@ void FeatureSpace::project_generated(double* prop, int size, std::vector<node_pt
             worst_score_ind = std::max_element(scores_sel.begin(), scores_sel.end()) - scores_sel.begin();
             for(int sc = 0; sc < scores_sel_private.size(); ++sc)
             {
-                if(((phi_sel.size() < _n_sis_select) || (scores_sel_private[sc] < scores_sel[worst_score_ind])) && _is_valid_feat_list(phi_sel_private[sc]->value_ptr(0), _n_samp, _cross_cor_max, phi_sel, scores_sel, scores_sel_private[sc]))
+                if(((phi_sel.size() < _n_sis_select) || (scores_sel_private[sc] < scores_sel[worst_score_ind])) && _is_valid_feat_list(phi_sel_private[sc]->value_ptr(), _n_samp, _cross_cor_max, phi_sel, scores_sel, scores_sel_private[sc]))
                 {
 
                     if(phi_sel.size() == _n_sis_select)
@@ -592,7 +588,6 @@ void FeatureSpace::sis(std::vector<double>& prop)
     {
         phi_sel.resize(cur_feat_local);
         scores_sel.resize(cur_feat_local);
-        node_ptr test_feat = phi_sel[0];
 
         project_generated(prop.data(), prop.size(), phi_sel, scores_sel);
 
@@ -640,7 +635,9 @@ void FeatureSpace::sis(std::vector<double>& prop)
                 {
                     sent_phi[inds[ii]]->set_value();
                     out_file_stream << std::setw(14) <<std::left << cur_feat + cur_feat_local << sent_phi[inds[ii]]->postfix_expr() << std::endl;
-                    sum_file_stream << std::setw(14) <<std::left << cur_feat + cur_feat_local << std::setw(24) << std::setprecision(18) << std::left << -1 * sent_scores[inds[ii]] << sent_phi[inds[ii]]->expr() << std::endl;
+//                    sum_file_stream << std::setw(14) <<std::left << cur_feat + cur_feat_local << std::setw(24) << std::setprecision(18) << std::left << -1 * sent_scores[inds[ii]] << '\t' << sent_phi[inds[ii]]->value_ptr()[0] << '\t' << sent_phi[inds[ii]]->feat(0)->value_ptr()[0] << '\t' << sent_phi[inds[ii]]->feat(0)->selected() << '\t' << sent_phi[inds[ii]]->feat(0)->arr_ind() << '\t' << node_value_arrs::N_STORE_FEATURES << sent_phi[inds[ii]]->feat(0)->expr() << '\t' << sent_phi[inds[ii]]->expr() << std::endl;
+
+                    sum_file_stream << std::setw(14) <<std::left << cur_feat + cur_feat_local << std::setw(18) << std::setprecision(12) << std::left << -1 * sent_scores[inds[ii]] << sent_phi[inds[ii]]->expr() << std::endl;
 
                     _phi_selected.push_back(sent_phi[inds[ii]]);
 
diff --git a/src/feature_creation/node/FeatureNode.hpp b/src/feature_creation/node/FeatureNode.hpp
index fed7e9f7c82643ed272809f824a9bf43c3e475cc..9002b0cac90b5a86ddfb200fc0d9a8e6f4d51bc6 100644
--- a/src/feature_creation/node/FeatureNode.hpp
+++ b/src/feature_creation/node/FeatureNode.hpp
@@ -167,7 +167,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    inline void set_value(int offset=0){std::copy_n(_value.data(), _n_samp, value_ptr());}
+    inline void set_value(int offset=0, bool for_comp=false){std::copy_n(_value.data(), _n_samp, value_ptr());}
 
     // DocString: feat_node_set_test_value
     /**
@@ -175,7 +175,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    inline void set_test_value(int offset=0){if(!_selected) std::copy_n(_test_value.data(), _n_test_samp, test_value_ptr());}
+    inline void set_test_value(int offset=0, bool for_comp=false){if(!_selected) std::copy_n(_test_value.data(), _n_test_samp, test_value_ptr());}
 
     // DocString: feat_node_is_nan
     /**
@@ -199,14 +199,14 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    inline double* value_ptr(int offset=0){return _selected ? node_value_arrs::get_d_matrix_ptr(_d_mat_ind) : node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, 0, offset);}
+    inline double* value_ptr(int offset=0, bool for_comp=false){return _selected ? node_value_arrs::get_d_matrix_ptr(_d_mat_ind) : node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, 0, offset, for_comp);}
 
     /**
      * @brief The pointer to where the feature's test data is stored
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    inline double* test_value_ptr(int offset=0){return node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, 0, offset);}
+    inline double* test_value_ptr(int offset=0, bool for_comp=false){return node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, 0, offset, for_comp);}
 
 
     // DocString: feat_node_rung
diff --git a/src/feature_creation/node/ModelNode.hpp b/src/feature_creation/node/ModelNode.hpp
index 4c65f29fb2ec119044d084628c99d0893eae727d..cc5ac624a2049b81849979d579bbcf9ec23f9ff5 100644
--- a/src/feature_creation/node/ModelNode.hpp
+++ b/src/feature_creation/node/ModelNode.hpp
@@ -116,7 +116,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    inline void set_value(int offset=0){return;}
+    inline void set_value(int offset=0, bool for_comp=false){return;}
 
     // DocString: model_node_set_test_value
     /**
@@ -124,7 +124,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    inline void set_test_value(int offset=0){return;}
+    inline void set_test_value(int offset=0, bool for_comp=false){return;}
 
     // DocString: model_node_is_nan
     /**
@@ -148,14 +148,14 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    inline double* value_ptr(int offset=0){return _value.data();}
+    inline double* value_ptr(int offset=0, bool for_comp=false){return _value.data();}
 
     /**
      * @brief The pointer to where the feature's test data is stored
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    inline double* test_value_ptr(int offset=0){return _test_value.data();}
+    inline double* test_value_ptr(int offset=0, bool for_comp=false){return _test_value.data();}
 
     // DocString: model_node_rung
     /**
diff --git a/src/feature_creation/node/Node.hpp b/src/feature_creation/node/Node.hpp
index 49f746883fd9be7c7c44d86a7c10a00718ac7532..3253e22d101ac7212998af9d5894ab4706c156ee 100644
--- a/src/feature_creation/node/Node.hpp
+++ b/src/feature_creation/node/Node.hpp
@@ -225,14 +225,14 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    virtual void set_value(int offset=0) = 0;
+    virtual void set_value(int offset=0, bool for_comp=false) = 0;
 
     /**
      * @brief The pointer to where the feature's training data is stored
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    virtual double* value_ptr(int offset=0) = 0;
+    virtual double* value_ptr(int offset=0, bool for_comp=false) = 0;
 
     // DocString: node_set_test_value
     /**
@@ -240,14 +240,14 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    virtual void set_test_value(int offset=0) = 0;
+    virtual void set_test_value(int offset=0, bool for_comp=false) = 0;
 
     /**
      * @brief The pointer to where the feature's test data is stored
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    virtual double* test_value_ptr(int offset=0) = 0;
+    virtual double* test_value_ptr(int offset=0, bool for_comp=false) = 0;
 
     // DocString: node_is_nan
     /**
diff --git a/src/feature_creation/node/operator_nodes/OperatorNode.hpp b/src/feature_creation/node/operator_nodes/OperatorNode.hpp
index a3c80b1bee7272b906f2ce297b481167dab70e28..4683df532dbcc546317a5792e7a55b1339ea0931 100644
--- a/src/feature_creation/node/operator_nodes/OperatorNode.hpp
+++ b/src/feature_creation/node/operator_nodes/OperatorNode.hpp
@@ -143,7 +143,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    virtual void set_value(int offset=0) = 0;
+    virtual void set_value(int offset=0, bool for_comp=false) = 0;
 
     // DocString: op_node_set_test_value
     /**
@@ -151,7 +151,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    virtual void set_test_value(int offset=0) = 0;
+    virtual void set_test_value(int offset=0, bool for_comp=false) = 0;
 
     /**
      * @brief Get the pointer to the feature's training data
@@ -161,17 +161,17 @@ public:
      *
      * @return pointer to the feature's training value
      */
-    double* value_ptr(int offset=0)
+    double* value_ptr(int offset=0, bool for_comp=false)
     {
         if(_selected)
             return node_value_arrs::get_d_matrix_ptr(_d_mat_ind);
 
-        if((rung() > node_value_arrs::N_RUNGS_STORED) && (node_value_arrs::temp_storage_reg(_arr_ind, rung(), offset) != _feat_ind))
+        if((rung() > node_value_arrs::N_RUNGS_STORED) && (node_value_arrs::temp_storage_reg(_arr_ind, rung(), offset, for_comp) != _feat_ind))
         {
-            set_value(offset);
+            set_value(offset, for_comp);
         }
 
-        return node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset);
+        return node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp);
     }
 
     /**
@@ -182,14 +182,14 @@ public:
      *
      * @return pointer to the feature's test values
      */
-    double* test_value_ptr(int offset=0)
+    double* test_value_ptr(int offset=0, bool for_comp=false)
     {
-        if((rung() > node_value_arrs::N_RUNGS_STORED) && (node_value_arrs::temp_storage_test_reg(_arr_ind, rung(), offset) != _feat_ind))
+        if((rung() > node_value_arrs::N_RUNGS_STORED) && (node_value_arrs::temp_storage_test_reg(_arr_ind, rung(), offset, for_comp) != _feat_ind))
         {
-            set_test_value(offset);
+            set_test_value(offset, for_comp);
         }
 
-        return node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset);
+        return node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp);
     }
 
     // DocString: op_node_is_nan
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs/absolute_value.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs/absolute_value.cpp
index 8807bb2b38567050846e00dc29208fa0ced00104..1d11bea55e722e2437f76cf180e5b2397baab114 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs/absolute_value.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs/absolute_value.cpp
@@ -63,15 +63,15 @@ void AbsNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     expected_abs_tot += std::abs(fact);
 }
 
-void AbsNode::set_value(int offset)
+void AbsNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::abs(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::abs(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::abs(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::abs(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void AbsNode::set_test_value(int offset)
+void AbsNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::abs(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::abs(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs/absolute_value.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs/absolute_value.hpp
index 25f3cd9dcc9e3df0265a0ddfae3917fdaeb4462a..901d5c1d31be9544526db29f15045f667eb0aeb7 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs/absolute_value.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs/absolute_value.hpp
@@ -86,7 +86,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: abs_node_set_test_value
     /**
@@ -94,7 +94,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: abs_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs_diff/absolute_difference.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs_diff/absolute_difference.cpp
index 7ec232c96023c58e0d990f2d3e67822d955be527..df01c20cf10ce1f7c8202684f6b5d1e35758ec84 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs_diff/absolute_difference.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs_diff/absolute_difference.cpp
@@ -84,15 +84,15 @@ void AbsDiffNode::update_div_mult_leaves(std::map<std::string, double>& div_mult
     expected_abs_tot += std::abs(fact);
 }
 
-void AbsDiffNode::set_value(int offset)
+void AbsDiffNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::abs_diff(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::abs_diff(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::abs_diff(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::abs_diff(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void AbsDiffNode::set_test_value(int offset)
+void AbsDiffNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::abs_diff(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), _feats[1]->test_value_ptr(2 * offset + 1), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::abs_diff(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), _feats[1]->test_value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs_diff/absolute_difference.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs_diff/absolute_difference.hpp
index 330c7512b5896961e72797919bdfe942b2dd24f7..0d6624b1cee6918ad68bfd6f008853653ee19a04 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs_diff/absolute_difference.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/abs_diff/absolute_difference.hpp
@@ -89,7 +89,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: abs_diff_node_set_test_value
     /**
@@ -97,7 +97,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: abs_diff_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add/add.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add/add.cpp
index dca9855669821a182a820c4969d717d043232db2..6df3da44a86961fd583113d2947a9fb7f30ab721 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add/add.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add/add.cpp
@@ -77,15 +77,15 @@ void AddNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     expected_abs_tot += std::abs(fact);
 }
 
-void AddNode::set_value(int offset)
+void AddNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::add(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::add(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::add(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::add(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void AddNode::set_test_value(int offset)
+void AddNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::add(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), _feats[1]->test_value_ptr(2 * offset + 1), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::add(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), _feats[1]->test_value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add/add.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add/add.hpp
index 3b3ca2faa5d8d4ab68bbc7aff020c663dbea9278..ee67f4421f50436b7fa797ce66a7b624b4a1e4ea 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add/add.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/add/add.hpp
@@ -86,7 +86,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: add_node_set_test_value
     /**
@@ -94,7 +94,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: add_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cb/cube.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cb/cube.cpp
index 1a849f02fdf0341082cbea93c36c4febcecb491e..f3088c0613de65dee8d252e2c42f6523b3a2a9e0 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cb/cube.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cb/cube.cpp
@@ -53,14 +53,14 @@ 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::set_value(int offset)
+void CbNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::cb(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::cb(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::cb(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::cb(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
-void CbNode::set_test_value(int offset)
+void CbNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::cb(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::cb(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cb/cube.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cb/cube.hpp
index ecb49fa79cea512a1c251f617665075ffc989f35..7fd9b2c7c62c0017b50967cd3a946fbebf3b3231 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cb/cube.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cb/cube.hpp
@@ -84,7 +84,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: cb_node_set_test_value
     /**
@@ -92,7 +92,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: cb_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cbrt/cube_root.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cbrt/cube_root.cpp
index 49a05fe909dd83c090a21ff93b4df3427dd9f922..9b156f38daf2e5048071dbd9314a2def9cef7fa8 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cbrt/cube_root.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cbrt/cube_root.cpp
@@ -52,15 +52,15 @@ 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::set_value(int offset)
+void CbrtNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::cbrt(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::cbrt(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::cbrt(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::cbrt(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void CbrtNode::set_test_value(int offset)
+void CbrtNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::cbrt(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::cbrt(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cbrt/cube_root.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cbrt/cube_root.hpp
index e0284b3ca8bd5848f832e705137c1e5e4806d7d1..516ab7a51b7fbf77812ef2cd4f46c31c0767e104 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cbrt/cube_root.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cbrt/cube_root.hpp
@@ -84,7 +84,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: cbrt_node_set_test_value
     /**
@@ -92,7 +92,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: cbrt_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos/cos.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos/cos.cpp
index 0cac2682551421e819dfebfa7064877b5458650c..3abacad261046077603bb13f3d19d492280f75ff 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos/cos.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos/cos.cpp
@@ -58,15 +58,15 @@ void CosNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     expected_abs_tot += std::abs(fact);
 }
 
-void CosNode::set_value(int offset)
+void CosNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::cos(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::cos(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::cos(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::cos(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void CosNode::set_test_value(int offset)
+void CosNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::cos(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::cos(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos/cos.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos/cos.hpp
index 523cb11b04e15f761ac9236eb002168266cb7130..96372f99ce39cbcbdd7617acc7d17f39db2e7494 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos/cos.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/cos/cos.hpp
@@ -84,7 +84,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: cos_node_set_test_value
     /**
@@ -92,7 +92,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: cos_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/div/divide.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/div/divide.cpp
index f4a0e005ec6643c17be8fb2838dadbb88b52d2d7..fe458088c60e87786df091692fcfe1b19895c1c4 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/div/divide.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/div/divide.cpp
@@ -78,15 +78,15 @@ void DivNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     _feats[1]->update_div_mult_leaves(div_mult_leaves, -1.0*fact, expected_abs_tot);
 }
 
-void DivNode::set_value(int offset)
+void DivNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::div(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::div(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::div(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::div(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void DivNode::set_test_value(int offset)
+void DivNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::div(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), _feats[1]->test_value_ptr(2 * offset + 1), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::div(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), _feats[1]->test_value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/div/divide.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/div/divide.hpp
index e5ef547b66ba04ab11a0de516b8136684ac363b6..870051e9029fa20a4d6aa85dc5f7baad97edb938 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/div/divide.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/div/divide.hpp
@@ -86,7 +86,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: div_node_set_test_value
     /**
@@ -94,7 +94,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: div_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exp/exponential.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exp/exponential.cpp
index d844a60fe62adeafcce09794af9ccee5ce719649..720b32a49b86de560ef665f080251fa70005cfef 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exp/exponential.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exp/exponential.cpp
@@ -60,15 +60,15 @@ void ExpNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     expected_abs_tot += std::abs(fact);
 }
 
-void ExpNode::set_value(int offset)
+void ExpNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::exp(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::exp(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::exp(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::exp(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void ExpNode::set_test_value(int offset)
+void ExpNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::exp(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::exp(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exp/exponential.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exp/exponential.hpp
index 2d17a56735a21073ee735ff624796ab6ca01d37a..1c87dc886d7325676179a16843c56e74c820afd6 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exp/exponential.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/exp/exponential.hpp
@@ -84,7 +84,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: exp_node_set_test_value
     /**
@@ -92,7 +92,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: exp_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inv/inverse.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inv/inverse.cpp
index aed5c35ec428269dc275f57e97147351a8472b4a..f9909ddebdc8b51bae581addb39a2e9d095295df 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inv/inverse.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inv/inverse.cpp
@@ -53,15 +53,15 @@ 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::set_value(int offset)
+void InvNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::inv(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::inv(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::inv(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::inv(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void InvNode::set_test_value(int offset)
+void InvNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::inv(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::inv(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inv/inverse.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inv/inverse.hpp
index 05c455e4b2ee0d668b4db1510546e11830ee81ac..ad32f06ca3b89d871960e790f8671a6d36722aaa 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inv/inverse.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/inv/inverse.hpp
@@ -72,7 +72,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: inv_node_set_test_value
     /**
@@ -80,7 +80,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: inv_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log/log.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log/log.cpp
index 0343bbd2243e3da70051aee4cf7eae73c839c743..3f05fa454b9e88cdf7c70e78b6b6d4387e22da96 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log/log.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log/log.cpp
@@ -61,15 +61,15 @@ void LogNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     expected_abs_tot += std::abs(fact);
 }
 
-void LogNode::set_value(int offset)
+void LogNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::log(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::log(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::log(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::log(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void LogNode::set_test_value(int offset)
+void LogNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::log(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::log(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log/log.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log/log.hpp
index d21bc141d00704ec769197789bca87b049d3b37e..35d2598fdd6f8634007e2154162a4fd9b0bc0aa3 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log/log.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/log/log.hpp
@@ -84,7 +84,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: log_node_set_test_value
     /**
@@ -92,7 +92,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: log_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/mult/multiply.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/mult/multiply.cpp
index 2641379ba3fc7231cdfb63f2ae90d1def6dfc8a0..5b0a88ba049c13aff3db289b04a88c66b617e924 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/mult/multiply.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/mult/multiply.cpp
@@ -78,15 +78,15 @@ void MultNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_le
     _feats[1]->update_div_mult_leaves(div_mult_leaves, fact, expected_abs_tot);
 }
 
-void MultNode::set_value(int offset)
+void MultNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::mult(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::mult(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::mult(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::mult(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void MultNode::set_test_value(int offset)
+void MultNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::mult(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), _feats[1]->test_value_ptr(2 * offset + 1), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::mult(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), _feats[1]->test_value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/mult/multiply.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/mult/multiply.hpp
index b6c73ff72365365bf407d4c92eac457961f928c6..1110bb2fdebefd79e9192b97248dcacf6ade8cfa 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/mult/multiply.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/mult/multiply.hpp
@@ -87,7 +87,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: mult_node_set_test_value
     /**
@@ -95,7 +95,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: mult_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/neg_exp/negative_exponential.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/neg_exp/negative_exponential.cpp
index 1ba3d169d80be38e2cbd85536d8dae33cea9bb2e..ca1fa72094c61aa27ca6da949f63a60438d71d73 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/neg_exp/negative_exponential.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/neg_exp/negative_exponential.cpp
@@ -61,15 +61,15 @@ void NegExpNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_
     expected_abs_tot += std::abs(fact);
 }
 
-void NegExpNode::set_value(int offset)
+void NegExpNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::neg_exp(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::neg_exp(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::neg_exp(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::neg_exp(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void NegExpNode::set_test_value(int offset)
+void NegExpNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::neg_exp(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::neg_exp(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/neg_exp/negative_exponential.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/neg_exp/negative_exponential.hpp
index 24880adffa8ff1e0f87867df5ded08aa14fec59c..3de95bebe5ef88e21add244f6e552b9420ae0076 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/neg_exp/negative_exponential.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/neg_exp/negative_exponential.hpp
@@ -85,7 +85,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: neg_exp_node_set_test_value
     /**
@@ -93,7 +93,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: neg_exp_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin/sin.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin/sin.cpp
index 1214263b275ba698bd2ac572f7893fed9b76995f..5fce04873fc789095d758f66d6924c2107d4c365 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin/sin.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin/sin.cpp
@@ -58,15 +58,15 @@ void SinNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     expected_abs_tot += std::abs(fact);
 }
 
-void SinNode::set_value(int offset)
+void SinNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::sin(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::sin(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::sin(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sin(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void SinNode::set_test_value(int offset)
+void SinNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::sin(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sin(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin/sin.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin/sin.hpp
index 3ad4f3de86d584db3462b28111f97483cbd9b108..500d785431bfb32aa3802673d46cf5fc13e35100 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin/sin.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sin/sin.hpp
@@ -85,7 +85,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: sin_node_set_test_value
     /**
@@ -93,7 +93,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: sin_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sp/sixth_power.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sp/sixth_power.cpp
index 3633e45675eb37407902444ce8236ec1693d5eed..000a5c993b7f4e5583185b808688ea2275052dbd 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sp/sixth_power.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sp/sixth_power.cpp
@@ -49,15 +49,15 @@ 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::set_value(int offset)
+void SixPowNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::sixth_pow(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::sixth_pow(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::sixth_pow(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sixth_pow(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void SixPowNode::set_test_value(int offset)
+void SixPowNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::sixth_pow(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sixth_pow(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sp/sixth_power.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sp/sixth_power.hpp
index 524bc61d10591aefc426cc1fcacbc40fa1e3fce7..180d11d13b042af629fcd20a78ec02833e199689 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sp/sixth_power.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sp/sixth_power.hpp
@@ -85,7 +85,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: six_pow_node_set_test_value
     /**
@@ -93,7 +93,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: six_pow_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sq/square.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sq/square.cpp
index 62d6ae37596dd4645684f206632ec3d30eb3df22..e4764ebaa8d60ba7fb7fb861854075b936ab8916 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sq/square.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sq/square.cpp
@@ -50,15 +50,15 @@ 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::set_value(int offset)
+void SqNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::sq(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::sq(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::sq(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sq(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void SqNode::set_test_value(int offset)
+void SqNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::sq(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sq(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sq/square.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sq/square.hpp
index b2e0097b479505504e43da06d9d8ce9d1ef015a4..81ed7b01419f595ea75e556f5615c6bf1a0944a2 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sq/square.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sq/square.hpp
@@ -84,7 +84,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: sq_node_set_test_value
     /**
@@ -92,7 +92,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: sq_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sqrt/square_root.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sqrt/square_root.cpp
index b549f3d5f17b1380dbd5acff534bd346d1953232..1ccdf3868a4108f9e3f7d684abe0671b0ed7b0c1 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sqrt/square_root.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sqrt/square_root.cpp
@@ -50,15 +50,15 @@ 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::set_value(int offset)
+void SqrtNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::sqrt(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::sqrt(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::sqrt(_n_samp, _feats[0]->value_ptr(2 * offset), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sqrt(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void SqrtNode::set_test_value(int offset)
+void SqrtNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::sqrt(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sqrt(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sqrt/square_root.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sqrt/square_root.hpp
index fb89ad2778fb60062dceb09aa133b2a6bf9c6904..f647601708add13c8cadaaedecfbfdb4191484d4 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sqrt/square_root.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sqrt/square_root.hpp
@@ -85,7 +85,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: sqrt_node_set_test_value
     /**
@@ -93,7 +93,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: sqrt_node_rung
     /**
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sub/subtract.cpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sub/subtract.cpp
index eae864d3778ab96abab73627b80ee4c84073b8da..7219048fbf683fd8f3f9cf21162083d4a6cef1a2 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sub/subtract.cpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sub/subtract.cpp
@@ -77,15 +77,15 @@ void SubNode::update_div_mult_leaves(std::map<std::string, double>& div_mult_lea
     expected_abs_tot += std::abs(fact);
 }
 
-void SubNode::set_value(int offset)
+void SubNode::set_value(int offset, bool for_comp)
 {
     if(_selected)
-        allowed_op_funcs::sub(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
+        allowed_op_funcs::sub(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_d_matrix_ptr(_d_mat_ind));
 
-    allowed_op_funcs::sub(_n_samp, _feats[0]->value_ptr(2 * offset), _feats[1]->value_ptr(2 * offset + 1), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sub(_n_samp, _feats[0]->value_ptr(2 * offset, for_comp), _feats[1]->value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
 
-void SubNode::set_test_value(int offset)
+void SubNode::set_test_value(int offset, bool for_comp)
 {
-    allowed_op_funcs::sub(_n_test_samp, _feats[0]->test_value_ptr(2 * offset), _feats[1]->test_value_ptr(2 * offset + 1), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset));
+    allowed_op_funcs::sub(_n_test_samp, _feats[0]->test_value_ptr(2 * offset, for_comp), _feats[1]->test_value_ptr(2 * offset + 1, for_comp), node_value_arrs::get_test_value_ptr(_arr_ind, _feat_ind, rung(), offset, for_comp));
 }
diff --git a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sub/subtract.hpp b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sub/subtract.hpp
index 54dc07b03ccd860d73b8b6a4ec0ab550c842c1dc..b20a1a512bc5f751712bad10a0d24a46088571a0 100644
--- a/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sub/subtract.hpp
+++ b/src/feature_creation/node/operator_nodes/allowed_operator_nodes/sub/subtract.hpp
@@ -87,7 +87,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_value(int offset=0);
+    void set_value(int offset=0, bool for_comp=false);
 
     // DocString: sub_node_set_test_value
     /**
@@ -95,7 +95,7 @@ public:
      *
      * @param offset(int) Key to determine which part of the temporary storage array to look into
      */
-    void set_test_value(int offset=0);
+    void set_test_value(int offset=0, bool for_comp=false);
 
     // DocString: sub_node_rung
     /**
diff --git a/src/feature_creation/node/value_storage/nodes_value_containers.cpp b/src/feature_creation/node/value_storage/nodes_value_containers.cpp
index 69eb30249e360034d54e7406409cd6fab5d4a397..da55a7d8014fa478ea7a4595818cc459a5ae3b54 100644
--- a/src/feature_creation/node/value_storage/nodes_value_containers.cpp
+++ b/src/feature_creation/node/value_storage/nodes_value_containers.cpp
@@ -40,7 +40,7 @@ void node_value_arrs::initialize_values_arr(int n_samples, int n_samples_test, i
     N_RUNGS_STORED = 0;
     N_STORE_FEATURES = n_primary_feat;
     MAX_RUNG = max_rung;
-    N_OP_SLOTS = static_cast<int>(std::pow(2, max_rung)) - 1;
+    N_OP_SLOTS = 2 * (static_cast<int>(std::pow(2, max_rung)) - 1);
 
     VALUES_ARR = std::vector<double>(N_STORE_FEATURES * N_SAMPLES);
     TEST_VALUES_ARR = std::vector<double>(N_STORE_FEATURES * N_SAMPLES_TEST);
@@ -98,7 +98,7 @@ void node_value_arrs::resize_values_arr(int n_dims, unsigned long int n_feat, bo
 
     if(use_temp)
     {
-        N_OP_SLOTS = static_cast<int>(std::pow(2, MAX_RUNG - N_RUNGS_STORED)) - 1;
+        N_OP_SLOTS = 2 * (static_cast<int>(std::pow(2, MAX_RUNG - N_RUNGS_STORED)) - 1);
 
         TEMP_STORAGE_ARR.resize(MAX_N_THREADS * (N_OP_SLOTS * N_STORE_FEATURES + 1) * N_SAMPLES);
         TEMP_STORAGE_ARR.shrink_to_fit();
@@ -122,22 +122,22 @@ void node_value_arrs::resize_values_arr(int n_dims, unsigned long int n_feat, bo
     }
 }
 
-double* node_value_arrs::get_value_ptr(unsigned long int arr_ind, unsigned long int feat_ind, int rung, int offset)
+double* node_value_arrs::get_value_ptr(unsigned long int arr_ind, unsigned long int feat_ind, int rung, int offset, bool for_comp)
 {
-    if(arr_ind < N_STORE_FEATURES)
+    if(rung <= N_RUNGS_STORED)
         return  access_value_arr(arr_ind);
-    int op_slot = get_op_slot(rung, offset);
 
+    int op_slot = get_op_slot(rung, offset, for_comp);
     temp_storage_reg(arr_ind, op_slot) = feat_ind;
     return access_temp_storage((arr_ind % N_STORE_FEATURES) + (op_slot % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1));
 }
 
-double* node_value_arrs::get_test_value_ptr(unsigned long int arr_ind, unsigned long int feat_ind, int rung, int offset)
+double* node_value_arrs::get_test_value_ptr(unsigned long int arr_ind, unsigned long int feat_ind, int rung, int offset, bool for_comp)
 {
-    if(arr_ind < N_STORE_FEATURES)
+    if(rung <= N_RUNGS_STORED)
         return  access_test_value_arr(arr_ind);
-    int op_slot = get_op_slot(rung, offset);
 
+    int op_slot = get_op_slot(rung, offset, for_comp);
     temp_storage_test_reg(arr_ind, op_slot) = feat_ind;
     return access_temp_storage_test((arr_ind % N_STORE_FEATURES) + (op_slot % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1));
 }
diff --git a/src/feature_creation/node/value_storage/nodes_value_containers.hpp b/src/feature_creation/node/value_storage/nodes_value_containers.hpp
index 4a2f088113559a059277c3a4dbc1161173dc139e..c00d37cd9cff2b97899cc486e7acdc691ddc0a43 100644
--- a/src/feature_creation/node/value_storage/nodes_value_containers.hpp
+++ b/src/feature_creation/node/value_storage/nodes_value_containers.hpp
@@ -15,6 +15,7 @@
 #include <memory>
 #include <numeric>
 #include <vector>
+#include <iostream>
 
 #include <omp.h>
 
@@ -126,7 +127,7 @@ namespace node_value_arrs
      *
      * @return The operator slot to use
      */
-    inline int get_op_slot(int rung, int offset){return std::abs(N_OP_SLOTS - static_cast<int>(std::pow(2, MAX_RUNG - rung)) - offset);}
+    inline int get_op_slot(int rung, int offset, bool for_comp){return std::abs(N_OP_SLOTS / (1 + !for_comp) - static_cast<int>(std::pow(2, MAX_RUNG - rung)) - offset);}
 
     /**
      * @brief Get a reference slot/feature register of the training data
@@ -136,7 +137,7 @@ namespace node_value_arrs
      *
      * @return The register element for a given feature index and op_slot
      */
-    inline int& temp_storage_reg(unsigned long int ind, int op_slot = 0){return TEMP_STORAGE_REG[(ind % N_STORE_FEATURES) + (op_slot % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1)];}
+    inline int& temp_storage_reg(unsigned long int ind, int op_slot=0){return TEMP_STORAGE_REG[(ind % N_STORE_FEATURES) + (op_slot % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1)];}
 
     /**
      * @brief Get a reference slot/feature register of the test data
@@ -146,7 +147,7 @@ namespace node_value_arrs
      *
      * @return The register element for a given feature index and op_slot
      */
-    inline int& temp_storage_test_reg(unsigned long int ind, int op_slot = 0){return TEMP_STORAGE_TEST_REG[(ind % N_STORE_FEATURES) + (op_slot % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1)];}
+    inline int& temp_storage_test_reg(unsigned long int ind, int op_slot=0){return TEMP_STORAGE_TEST_REG[(ind % N_STORE_FEATURES) + (op_slot % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1)];}
 
     /**
      * @brief Get a reference slot/feature register of the training data
@@ -157,7 +158,7 @@ namespace node_value_arrs
      *
      * @return The register element for a given feature index and offset
      */
-    inline int& temp_storage_reg(unsigned long int ind, int rung, int offset){return TEMP_STORAGE_REG[(ind % N_STORE_FEATURES) + (get_op_slot(rung, offset) % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1)];}
+    inline int& temp_storage_reg(unsigned long int ind, int rung, int offset, bool for_comp){return TEMP_STORAGE_REG[(ind % N_STORE_FEATURES) + (get_op_slot(rung, offset, for_comp) % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1)];}
 
     /**
      * @brief Get a reference slot/feature register of the test data
@@ -168,7 +169,7 @@ namespace node_value_arrs
      *
      * @return The register element for a given feature index and offset
      */
-    inline int& temp_storage_test_reg(unsigned long int ind, int rung, int offset){return TEMP_STORAGE_TEST_REG[(ind % N_STORE_FEATURES) + (get_op_slot(rung, offset) % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1)];}
+    inline int& temp_storage_test_reg(unsigned long int ind, int rung, int offset, bool for_comp){return TEMP_STORAGE_TEST_REG[(ind % N_STORE_FEATURES) + (get_op_slot(rung, offset, for_comp) % N_OP_SLOTS) * N_STORE_FEATURES + omp_get_thread_num() * (N_STORE_FEATURES * N_OP_SLOTS + 1)];}
 
     /**
      * @brief Access element of the permanent training data storage array
@@ -177,7 +178,7 @@ namespace node_value_arrs
      *
      * @return pointer to the Node's training data
      */
-    inline double* access_value_arr(unsigned long int feature_ind){return VALUES_ARR.data() + feature_ind*N_SAMPLES;}
+    inline double* access_value_arr(unsigned long int feature_ind){return &VALUES_ARR[feature_ind*N_SAMPLES];}
 
     /**
      * @brief Access element of the permanent test data storage array
@@ -186,7 +187,7 @@ namespace node_value_arrs
      *
      * @return pointer to the Node's test data
      */
-    inline double* access_test_value_arr(unsigned long int feature_ind){return TEST_VALUES_ARR.data() + feature_ind*N_SAMPLES_TEST;}
+    inline double* access_test_value_arr(unsigned long int feature_ind){return &TEST_VALUES_ARR[feature_ind*N_SAMPLES_TEST];}
 
     /**
      * @brief Access element of temporary storage array for the training data
@@ -195,7 +196,7 @@ namespace node_value_arrs
      *
      * @return pointer to the data stored in the specified slot
      */
-    inline double* access_temp_storage(unsigned long int slot){return TEMP_STORAGE_ARR.data() + slot*N_SAMPLES;}
+    inline double* access_temp_storage(unsigned long int slot){return &TEMP_STORAGE_ARR[slot*N_SAMPLES];}
 
     /**
      * @brief Access element of temporary storage array for the test data
@@ -204,7 +205,7 @@ namespace node_value_arrs
      *
      * @return pointer to the data stored in the specified slot
      */
-    inline double* access_temp_storage_test(unsigned long int slot){return TEMP_STORAGE_TEST_ARR.data() + slot*N_SAMPLES_TEST;}
+    inline double* access_temp_storage_test(unsigned long int slot){return &TEMP_STORAGE_TEST_ARR[slot*N_SAMPLES_TEST];}
 
     /**
      * @brief Get a Node's value_ptr
@@ -215,7 +216,7 @@ namespace node_value_arrs
      *
      * @return The value pointer
      */
-    double* get_value_ptr(unsigned long int arr_ind, unsigned long int feat_ind, int rung = 0, int offset = 0);
+    double* get_value_ptr(unsigned long int arr_ind, unsigned long int feat_ind, int rung=0, int offset=0, bool for_comp=false);
 
     /**
      * @brief Get a Node's test_value_ptr
@@ -226,7 +227,7 @@ namespace node_value_arrs
      *
      * @return The value pointer
      */
-    double* get_test_value_ptr(unsigned long int arr_ind, unsigned long int feat_ind, int rung = 0, int offset = 0);
+    double* get_test_value_ptr(unsigned long int arr_ind, unsigned long int feat_ind, int rung=0, int offset=0, bool for_comp=false);
 
     /**
      * @brief Get the pointer to a particular selected Node from sis
@@ -234,28 +235,25 @@ namespace node_value_arrs
      * @param ind Index of the data in the descriptor matrix
      * @return The pointer to the descriptor matrix's data
      */
-    inline double* get_d_matrix_ptr(int ind)
-    {
-        return D_MATRIX.data() + ind * N_SAMPLES;
-    }
+    inline double* get_d_matrix_ptr(int ind){return &D_MATRIX[ind * N_SAMPLES];}
 
     /**
      * @brief Flush the temporary storage register (training data)
      * @details Reset all slots in the register to -1
      */
-    inline void clear_temp_reg(){std::fill_n(TEMP_STORAGE_REG.data(), TEMP_STORAGE_REG.size(), -1);}
+    inline void clear_temp_reg(){std::fill_n(TEMP_STORAGE_REG.begin(), TEMP_STORAGE_REG.size(), -1);}
 
     /**
      * @brief Flush the temporary storage register (training data)
      * @details Reset all slots in the register to -1
      */
-    inline void clear_temp_reg_thread(){std::fill_n(TEMP_STORAGE_REG.data() + N_STORE_FEATURES * N_OP_SLOTS * omp_get_thread_num(), N_STORE_FEATURES * N_OP_SLOTS, -1);}
+    inline void clear_temp_reg_thread(){std::fill_n(TEMP_STORAGE_REG.begin() + (N_STORE_FEATURES * N_OP_SLOTS + 1) * omp_get_thread_num(), N_STORE_FEATURES * N_OP_SLOTS + 1, -1);}
 
     /**
      * @brief Flush the temporary storage register (test data)
      * @details Reset all slots in the register to -1
      */
-    inline void clear_temp_test_reg(){std::fill_n(TEMP_STORAGE_TEST_REG.data(), TEMP_STORAGE_TEST_REG.size(), -1);}
+    inline void clear_temp_test_reg(){std::fill_n(TEMP_STORAGE_TEST_REG.begin(), TEMP_STORAGE_TEST_REG.size(), -1);}
 
 }
 
diff --git a/src/python/bindings_docstring_keyed.hpp b/src/python/bindings_docstring_keyed.hpp
index 1687d8074afc0f9a7ce48cdeedef87571709c914..4e3da800083d0f2bc11bb0cda94455d79b5f80f8 100644
--- a/src/python/bindings_docstring_keyed.hpp
+++ b/src/python/bindings_docstring_keyed.hpp
@@ -42,10 +42,10 @@ namespace sisso
                 inline Unit unit(){return this->get_override("unit")();}
                 inline std::vector<double> value(){return this->get_override("value")();}
                 inline std::vector<double> test_value(){return this->get_override("test_value")();}
-                inline void set_value(int offset=0){this->get_override("set_value")();}
-                inline double* value_ptr(int offset=0){return this->get_override("value_ptr")();}
-                inline void set_test_value(int offset=0){this->get_override("set_test_value")();}
-                inline double* test_value_ptr(int offset=0){return this->get_override("test_value_ptr")();}
+                inline void set_value(int offset=0, bool for_comp=false){this->get_override("set_value")();}
+                inline double* value_ptr(int offset=0, bool for_comp=false){return this->get_override("value_ptr")();}
+                inline void set_test_value(int offset=0, bool for_comp=false){this->get_override("set_test_value")();}
+                inline double* test_value_ptr(int offset=0, bool for_comp=false){return this->get_override("test_value_ptr")();}
                 inline bool is_nan(){return this->get_override("is_nan")();}
                 inline bool is_const(){return this->get_override("is_const")();}
                 inline NODE_TYPE type(){return this->get_override("type")();}
@@ -65,8 +65,8 @@ namespace sisso
             template<int N>
             struct OperatorNodeWrap : OperatorNode<N>, py::wrapper<OperatorNode<N>>
             {
-                inline void set_value(int offset=0){this->get_override("set_value")();}
-                inline void set_test_value(int offset=0){this->get_override("set_test_value")();}
+                inline void set_value(int offset=0, bool for_comp=false){this->get_override("set_value")();}
+                inline void set_test_value(int offset=0, bool for_comp=false){this->get_override("set_test_value")();}
                 inline NODE_TYPE type(){return this->get_override("type")();}
                 inline int rung(int cur_rung = 0){return this->get_override("rung")();}
                 inline std::string expr(){return this->get_override("expr")();}
diff --git a/src/utils/compare_features.cpp b/src/utils/compare_features.cpp
index 5eabd40018299566528ef5de0a4e65b312e1f909..2720345cbf269fc920163992e5bea4d08b26e4ec 100644
--- a/src/utils/compare_features.cpp
+++ b/src/utils/compare_features.cpp
@@ -1,5 +1,5 @@
 #include <utils/compare_features.hpp>
-
+#include <iostream>
 std::vector<double> comp_feats::RANK;
 std::vector<int> comp_feats::INDEX;
 
@@ -44,14 +44,15 @@ void comp_feats::set_is_valid_fxn(
 
 bool comp_feats::valid_feature_against_selected_max_corr_1_pearson(double* val_ptr, int n_samp, double cross_cor_max, std::vector<double>& scores_sel, double cur_score, int end_sel, int start_sel)
 {
-    double base_val = util_funcs::r(val_ptr, val_ptr, n_samp);
+    double mean = util_funcs::mean(val_ptr, n_samp);
+    double stand_dev = util_funcs::stand_dev(val_ptr, n_samp, mean);
+    double base_val = util_funcs::r(val_ptr, val_ptr, n_samp, mean, stand_dev, mean, stand_dev);
 
     for(int dd = start_sel; dd < end_sel; ++dd)
     {
-        if(abs(cur_score - scores_sel[dd]) > 1e-5)
+        if(std::abs(cur_score - scores_sel[dd]) > 1e-5)
             continue;
-
-        if((base_val - std::abs(util_funcs::r(node_value_arrs::get_d_matrix_ptr(dd), val_ptr, n_samp))) < 1e-9)
+        if((base_val - std::abs(util_funcs::r(val_ptr, node_value_arrs::get_d_matrix_ptr(dd), n_samp, mean, stand_dev))) < 5e-9)
             return false;
     }
     return true;
@@ -59,14 +60,16 @@ bool comp_feats::valid_feature_against_selected_max_corr_1_pearson(double* val_p
 
 bool comp_feats::valid_feature_against_selected_max_corr_1_feat_list_pearson(double* val_ptr, int n_samp, double cross_cor_max, std::vector<node_ptr>& selected, std::vector<double>& scores_sel, double cur_score)
 {
-    double base_val = util_funcs::r(val_ptr, val_ptr, n_samp);
+    double mean = util_funcs::mean(val_ptr, n_samp);
+    double stand_dev = util_funcs::stand_dev(val_ptr, n_samp, mean);
+    double base_val = util_funcs::r(val_ptr, val_ptr, n_samp, mean, stand_dev, mean, stand_dev);
 
     for(int ff = 0; ff < selected.size(); ++ff)
     {
-        if(abs(cur_score - scores_sel[ff]) > 1e-5)
+        if(std::abs(cur_score - scores_sel[ff]) > 1e-5)
             continue;
 
-        if((base_val - std::abs(util_funcs::r(selected[ff]->value_ptr(1), val_ptr, n_samp))) < 1e-9)
+        if((base_val - std::abs(util_funcs::r(val_ptr, selected[ff]->value_ptr(0, true), n_samp, mean, stand_dev))) < 5e-9)
             return false;
     }
     return true;
@@ -75,7 +78,10 @@ bool comp_feats::valid_feature_against_selected_max_corr_1_feat_list_pearson(dou
 
 bool comp_feats::valid_feature_against_selected_pearson(double* val_ptr, int n_samp, double cross_cor_max, std::vector<double>& scores_sel, double cur_score, int end_sel, int start_sel)
 {
-    double base_val = util_funcs::r(val_ptr, val_ptr, n_samp);
+    double mean = util_funcs::mean(val_ptr, n_samp);
+    double stand_dev = util_funcs::stand_dev(val_ptr, n_samp, mean);
+    double base_val = util_funcs::r(val_ptr, val_ptr, n_samp, mean, stand_dev, mean, stand_dev);
+
     volatile bool is_valid = true;
 
     #pragma omp parallel for schedule(dynamic)
@@ -84,7 +90,7 @@ bool comp_feats::valid_feature_against_selected_pearson(double* val_ptr, int n_s
         if(!is_valid)
             continue;
 
-        if((base_val - std::abs(util_funcs::r(node_value_arrs::get_d_matrix_ptr(dd), val_ptr, n_samp))) < (1.0 - cross_cor_max + 1e-10))
+        if((base_val - std::abs(util_funcs::r(val_ptr, node_value_arrs::get_d_matrix_ptr(dd), n_samp, mean, stand_dev))) < (1.0 - cross_cor_max + 5e-9))
             is_valid = false;
     }
     return is_valid;
@@ -92,11 +98,13 @@ bool comp_feats::valid_feature_against_selected_pearson(double* val_ptr, int n_s
 
 bool comp_feats::valid_feature_against_selected_feat_list_pearson(double* val_ptr, int n_samp, double cross_cor_max, std::vector<node_ptr>& selected, std::vector<double>& scores_sel, double cur_score)
 {
-    double base_val = util_funcs::r(val_ptr, val_ptr, n_samp);
+    double mean = util_funcs::mean(val_ptr, n_samp);
+    double stand_dev = util_funcs::stand_dev(val_ptr, n_samp, mean);
+    double base_val = util_funcs::r(val_ptr, val_ptr, n_samp, mean, stand_dev, mean, stand_dev);
 
     for(auto& feat : selected)
     {
-        if((base_val - std::abs(util_funcs::r(feat->value_ptr(1), val_ptr, n_samp))) < (1.0 - cross_cor_max + 1e-10))
+        if((base_val - std::abs(util_funcs::r(val_ptr, feat->value_ptr(0, true), n_samp, mean, stand_dev))) < (1.0 - cross_cor_max + 5e-9))
             return false;
     }
     return true;
@@ -108,12 +116,12 @@ bool comp_feats::valid_feature_against_selected_max_corr_1_spearman(double* val_
 
     for(int dd = start_sel; dd < end_sel; ++dd)
     {
-        if(abs(std::floor(cur_score) - std::floor(scores_sel[dd])) > 1e-5)
+        if(std::abs(std::floor(cur_score) - std::floor(scores_sel[dd])) > 1e-5)
             continue;
 
         // Rank the new variable and take the Pearson correlation of the rank variables (val_ptr rank still in &RANK[(omp_get_thread_num() * 4 + 2) * n_samp])
         util_funcs::rank(node_value_arrs::get_d_matrix_ptr(dd), &RANK[omp_get_thread_num() * 4 * n_samp], &INDEX[omp_get_thread_num() * 2 * n_samp], n_samp);
-        if((base_val - std::abs(util_funcs::r(&RANK[omp_get_thread_num() * 4 * n_samp], &RANK[(omp_get_thread_num() * 4 + 2) * n_samp], n_samp))) < 1e-9)
+        if((base_val - std::abs(util_funcs::r(&RANK[omp_get_thread_num() * 4 * n_samp], &RANK[(omp_get_thread_num() * 4 + 2) * n_samp], n_samp))) < 5e-9)
             return false;
     }
     return true;
@@ -125,12 +133,12 @@ bool comp_feats::valid_feature_against_selected_max_corr_1_feat_list_spearman(do
 
     for(int ff = 0; ff < selected.size(); ++ff)
     {
-        if(abs(std::floor(cur_score) - std::floor(scores_sel[ff])) > 1e-5)
+        if(std::abs(std::floor(cur_score) - std::floor(scores_sel[ff])) > 1e-5)
             continue;
 
         // Rank the new variable and take the Pearson correlation of the rank variables (val_ptr rank still in &RANK[(omp_get_thread_num() * 4 + 2) * n_samp])
-        util_funcs::rank(selected[ff]->value_ptr(1), &RANK[omp_get_thread_num() * 4 * n_samp], &INDEX[omp_get_thread_num() * 2 * n_samp], n_samp);
-        if((base_val - std::abs(util_funcs::r(&RANK[omp_get_thread_num() * 4 * n_samp], &RANK[(omp_get_thread_num() * 4 + 2) * n_samp], n_samp))) < 1e-9)
+        util_funcs::rank(selected[ff]->value_ptr(0, true), &RANK[omp_get_thread_num() * 4 * n_samp], &INDEX[omp_get_thread_num() * 2 * n_samp], n_samp);
+        if((base_val - std::abs(util_funcs::r(&RANK[omp_get_thread_num() * 4 * n_samp], &RANK[(omp_get_thread_num() * 4 + 2) * n_samp], n_samp))) < 5e-9)
             return false;
     }
     return true;
@@ -150,7 +158,7 @@ bool comp_feats::valid_feature_against_selected_spearman(double* val_ptr, int n_
 
         // Rank the new variable and take the Pearson correlation of the rank variables (val_ptr rank still in &RANK[(omp_get_thread_num() * 4 + 2) * n_samp])
         util_funcs::rank(node_value_arrs::get_d_matrix_ptr(dd), &RANK[omp_get_thread_num() * 4 * n_samp], &INDEX[omp_get_thread_num() * 2 * n_samp], n_samp);
-        if((base_val - std::abs(util_funcs::r(&RANK[omp_get_thread_num() * 4 * n_samp], &RANK[(omp_get_thread_num() * 4 + 2) * n_samp], n_samp))) < cross_cor_max + 1e-9)
+        if((base_val - std::abs(util_funcs::r(&RANK[omp_get_thread_num() * 4 * n_samp], &RANK[(omp_get_thread_num() * 4 + 2) * n_samp], n_samp))) < cross_cor_max + 5e-9)
             is_valid = false;
     }
     return is_valid;
@@ -163,8 +171,8 @@ bool comp_feats::valid_feature_against_selected_feat_list_spearman(double* val_p
     for(auto& feat : selected)
     {
         // Rank the new variable and take the Pearson correlation of the rank variables (val_ptr rank still in &RANK[(omp_get_thread_num() * 4 + 2) * n_samp])
-        util_funcs::rank(feat->value_ptr(1), &RANK[omp_get_thread_num() * 4 * n_samp], &INDEX[omp_get_thread_num() * 2 * n_samp], n_samp);
-        if((base_val - std::abs(util_funcs::r(&RANK[omp_get_thread_num() * 4 * n_samp], &RANK[(omp_get_thread_num() * 4 + 2) * n_samp], n_samp))) < cross_cor_max + 1e-9)
+        util_funcs::rank(feat->value_ptr(0, true), &RANK[omp_get_thread_num() * 4 * n_samp], &INDEX[omp_get_thread_num() * 2 * n_samp], n_samp);
+        if((base_val - std::abs(util_funcs::r(&RANK[omp_get_thread_num() * 4 * n_samp], &RANK[(omp_get_thread_num() * 4 + 2) * n_samp], n_samp))) < cross_cor_max + 5e-9)
             return false;
     }
     return true;
diff --git a/src/utils/math_funcs.cpp b/src/utils/math_funcs.cpp
index 2b1fd04169a3687673d06bbbf6268906831dcdcb..2d8699546962f0150de655d75569e6a5b2146111 100644
--- a/src/utils/math_funcs.cpp
+++ b/src/utils/math_funcs.cpp
@@ -24,6 +24,12 @@ double util_funcs::log_r2(double* a, double* b, double* log_a, int size)
     return r2(log_a, b, size);
 }
 
+double util_funcs::log_r2(double* a, double* b, double* log_a, int size, double mean_b, double std_b)
+{
+    std::transform(a, a + size, log_a, [](double aa){return std::log(aa);});
+    return r2(b, log_a, size, mean_b, std_b);
+}
+
 double util_funcs::r(double* a, double* b, std::vector<int>& sizes)
 {
     double result = 0.0;
@@ -37,6 +43,32 @@ double util_funcs::r(double* a, double* b, std::vector<int>& sizes)
     return std::sqrt(result / sizes.size());
 }
 
+double util_funcs::r(double* a, double* b, std::vector<int>& sizes, std::vector<double>& mean_a, std::vector<double>& std_a)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < sizes.size(); ++nt)
+    {
+        result += r2(a + pos, b + pos, sizes[nt], mean_a[nt], std_a[nt]);
+        pos += sizes[nt];
+    }
+
+    return std::sqrt(result / sizes.size());
+}
+
+double util_funcs::r(double* a, double* b, std::vector<int>& sizes, std::vector<double>& mean_a, std::vector<double>& std_a, std::vector<double>& mean_b, std::vector<double>& std_b)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < sizes.size(); ++nt)
+    {
+        result += r2(a + pos, b + pos, sizes[nt], mean_a[nt], std_a[nt], mean_b[nt], std_b[nt]);
+        pos += sizes[nt];
+    }
+
+    return std::sqrt(result / sizes.size());
+}
+
 double util_funcs::r2(double* a, double* b, std::vector<int>& sizes)
 {
     double result = 0.0;
@@ -50,6 +82,32 @@ double util_funcs::r2(double* a, double* b, std::vector<int>& sizes)
     return result / sizes.size();
 }
 
+double util_funcs::r2(double* a, double* b, std::vector<int>& sizes, std::vector<double>& mean_a, std::vector<double>& std_a)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < sizes.size(); ++nt)
+    {
+        result += r2(a + pos, b + pos, sizes[nt], mean_a[nt], std_a[nt]);
+        pos += sizes[nt];
+    }
+
+    return result / sizes.size();
+}
+
+double util_funcs::r2(double* a, double* b, std::vector<int>& sizes, std::vector<double>& mean_a, std::vector<double>& std_a, std::vector<double>& mean_b, std::vector<double>& std_b)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < sizes.size(); ++nt)
+    {
+        result += r2(a + pos, b + pos, sizes[nt], mean_a[nt], std_a[nt], mean_b[nt], std_b[nt]);
+        pos += sizes[nt];
+    }
+
+    return result / sizes.size();
+}
+
 double util_funcs::log_r2(double* a, double* b, double* log_a, std::vector<int>& sizes)
 {
     double result = 0.0;
@@ -63,6 +121,19 @@ double util_funcs::log_r2(double* a, double* b, double* log_a, std::vector<int>&
     return result / sizes.size();
 }
 
+double util_funcs::log_r2(double* a, double* b, double* log_a, std::vector<int>& sizes, std::vector<double>& mean_b, std::vector<double>& std_b)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < sizes.size(); ++nt)
+    {
+        result += log_r2(a + pos, b + pos, log_a, sizes[nt], mean_b[nt], std_b[nt]);
+        pos += sizes[nt];
+    }
+
+    return result / sizes.size();
+}
+
 double util_funcs::r(double* a, double* b, int* sz, int n_task)
 {
     double result = 0.0;
@@ -76,6 +147,32 @@ double util_funcs::r(double* a, double* b, int* sz, int n_task)
     return std::sqrt(result / static_cast<double>(n_task));
 }
 
+double util_funcs::r(double* a, double* b, int* sz, double* mean_a, double* std_a, int n_task)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < n_task; ++nt)
+    {
+        result += r2(a + pos, b + pos, sz[nt], mean_a[nt], std_a[nt]);
+        pos += sz[nt];
+    }
+
+    return std::sqrt(result / static_cast<double>(n_task));
+}
+
+double util_funcs::r(double* a, double* b, int* sz, double* mean_a, double* std_a, double* mean_b, double* std_b, int n_task)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < n_task; ++nt)
+    {
+        result += r2(a + pos, b + pos, sz[nt], mean_a[nt], std_a[nt], mean_b[nt], std_b[nt]);
+        pos += sz[nt];
+    }
+
+    return std::sqrt(result / static_cast<double>(n_task));
+}
+
 double util_funcs::r2(double* a, double* b, int* sz, int n_task)
 {
     double result = 0.0;
@@ -89,6 +186,32 @@ double util_funcs::r2(double* a, double* b, int* sz, int n_task)
     return result / static_cast<double>(n_task);
 }
 
+double util_funcs::r2(double* a, double* b, int* sz, double* mean_a, double* std_a, int n_task)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < n_task; ++nt)
+    {
+        result += r2(a + pos, b + pos, sz[nt], mean_a[nt], std_a[nt]);
+        pos += sz[nt];
+    }
+
+    return result / static_cast<double>(n_task);
+}
+
+double util_funcs::r2(double* a, double* b, int* sz, double* mean_a, double* std_a, double* mean_b, double* std_b, int n_task)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < n_task; ++nt)
+    {
+        result += r2(a + pos, b + pos, sz[nt], mean_a[nt], std_a[nt], mean_b[nt], std_b[nt]);
+        pos += sz[nt];
+    }
+
+    return result / static_cast<double>(n_task);
+}
+
 double util_funcs::log_r2(double* a, double* b, double* log_a, int* sz, int n_task)
 {
     double result = 0.0;
@@ -102,6 +225,19 @@ double util_funcs::log_r2(double* a, double* b, double* log_a, int* sz, int n_ta
     return result / static_cast<double>(n_task);
 }
 
+double util_funcs::log_r2(double* a, double* b, double* log_a, int* sz, double* mean_b, double* std_b, int n_task)
+{
+    double result = 0.0;
+    int pos = 0;
+    for(int nt = 0; nt < n_task; ++nt)
+    {
+        result += log_r2(a + pos, b + pos, log_a + pos, sz[nt], mean_b[nt], std_b[nt]);
+        pos += sz[nt];
+    }
+
+    return result / static_cast<double>(n_task);
+}
+
 std::vector<int> util_funcs::argsort(std::vector<double>& vec)
 {
     std::vector<int> index(vec.size());
diff --git a/src/utils/math_funcs.hpp b/src/utils/math_funcs.hpp
index 5df0d7e3dcc3b6fc4f6b5914ef72685f6e44cdbd..9965cb9781d6bfa9a9b0b0264aec32210abf5230 100644
--- a/src/utils/math_funcs.hpp
+++ b/src/utils/math_funcs.hpp
@@ -132,11 +132,46 @@ namespace util_funcs
      * @param a the pointer to the head of the first vector
      * @param b the pointer to the head of the second vector
      * @param size the size of the vector
-     * @return [description]
+     * @return The correlation coefficient between vector a and vector b
      */
     inline double r(double* a, double* b, int size)
     {
-        return std::inner_product(a, a + size, b, -1.0 * static_cast<double>(size) * mean(a, size) * mean(b, size)) / (static_cast<double>(size) * stand_dev(a, size) * stand_dev(b, size));
+        double mean_a = mean(a, size);
+        double mean_b = mean(b, size);
+        return std::inner_product(a, a + size, b, -1.0 * static_cast<double>(size) * mean_a * mean_b) / (static_cast<double>(size) * stand_dev(a, size, mean_a) * stand_dev(b, size, mean_b));
+    }
+
+    /**
+     * @brief The Pearson correlation for two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param size the size of the vector
+     * @param mean_a the mean of the a vector
+     * @param std_a the standard deviation of the a vector
+     * @return The correlation coefficient between vector a and vector b
+     */
+    inline double r(double* a, double* b, int size, double mean_a, double std_a)
+    {
+        double mean_b = mean(b, size);
+        return std::inner_product(a, a + size, b, -1.0 * static_cast<double>(size) * mean_a * mean_b) / (static_cast<double>(size) * std_a * stand_dev(b, size, mean_b));
+    }
+
+    /**
+     * @brief The Pearson correlation for two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param size the size of the vector
+     * @param mean_a the mean of the a vector
+     * @param std_a the standard deviation of the a vector
+     * @param mean_b the mean of the b vector
+     * @param std_b the standard deviation of the b vector
+     * @return The correlation coefficient between vector a and vector b
+     */
+    inline double r(double* a, double* b, int size, double mean_a, double std_a, double mean_b, double std_b)
+    {
+        return std::inner_product(a, a + size, b, -1.0 * static_cast<double>(size) * mean_a * mean_b) / (static_cast<double>(size) * std_a * std_b);
     }
 
     /**
@@ -149,6 +184,32 @@ namespace util_funcs
      */
     double r(double* a, double* b, std::vector<int>& sizes);
 
+    /**
+     * @brief The Pearson correlation for two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param sizes the sizes of the tasks to calculate the correlation on
+     * @param mean_a the mean of the a vector for each task
+     * @param std_a the standard deviation of the a vector for each task
+     * @return The average Pearson correlations
+     */
+    double r(double* a, double* b, std::vector<int>& sizes, std::vector<double>& mean_a, std::vector<double>& std_a);
+
+    /**
+     * @brief The Pearson correlation for two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param sizes the sizes of the tasks to calculate the correlation on
+     * @param mean_a the mean of the a vector for each task
+     * @param std_a the standard deviation of the a vector for each task
+     * @param mean_b the mean of the b vector for each task
+     * @param std_b the standard deviation of the b vector for each task
+     * @return The average Pearson correlations
+     */
+    double r(double* a, double* b, std::vector<int>& sizes, std::vector<double>& mean_a, std::vector<double>& std_a, std::vector<double>& mean_b, std::vector<double>& std_b);
+
     /**
      * @brief The Pearson correlation for two vectors
      *
@@ -161,16 +222,79 @@ namespace util_funcs
     double r(double* a, double* b, int* sz, int n_tasks);
 
     /**
-     * @brief Calculate the Coefficient of Determination between two vectors
+     * @brief The Pearson correlation for two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param sz the start of vector that describes the sizes of the tasks to calculate the correlation on
+     * @param mean_a the mean of the a vector for each task
+     * @param std_a the standard deviation of the a vector for each task
+     * @param n_tasks number of tasks to average over
+     * @return The average Pearson correlations
+     */
+    double r(double* a, double* b, int* sz, double* mean_a, double* std_a, int n_tasks);
+
+    /**
+     * @brief The Pearson correlation for two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param sz the start of vector that describes the sizes of the tasks to calculate the correlation on
+     * @param mean_a the mean of the a vector for each task
+     * @param std_a the standard deviation of the b vector for each task
+     * @param mean_b the mean of the a vector for each task
+     * @param std_b the standard deviation of the b vector for each task
+     * @param n_tasks number of tasks to average over
+     * @return The average Pearson correlations
+     */
+    double r(double* a, double* b, int* sz, double* mean_a, double* std_a, double* mean_b, double* std_b, int n_tasks);
+
+    /**
+     * @brief The Pearson correlation for two vectors
      *
      * @param a the pointer to the head of the first vector
      * @param b the pointer to the head of the second vector
      * @param size the size of the vector
-     * @return The Coefficient of Determination
+     * @return The coefficient of determination between vector a and vector b
      */
     inline double r2(double* a, double* b, int size)
     {
-        return std::pow(std::inner_product(a, a + size, b, -1.0 * static_cast<double>(size) * mean(a, size) * mean(b, size)) / (static_cast<double>(size) * stand_dev(a, size) * stand_dev(b, size)), 2.0);
+        double mean_a = mean(a, size);
+        double mean_b = mean(b, size);
+        return std::pow(std::inner_product(a, a + size, b, -1.0 * static_cast<double>(size) * mean_a * mean_b) / (static_cast<double>(size) * stand_dev(a, size, mean_a) * stand_dev(b, size, mean_b)), 2.0);
+    }
+
+    /**
+     * @brief The Pearson correlation for two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param size the size of the vector
+     * @param mean_a the mean of the a vector
+     * @param std_a the standard deviation of the a vector
+     * @return The coefficient of determination between vector a and vector b
+     */
+    inline double r2(double* a, double* b, int size, double mean_a, double std_a)
+    {
+        double mean_b = mean(b, size);
+        return std::pow(std::inner_product(a, a + size, b, -1.0 * static_cast<double>(size) * mean_a * mean_b) / (static_cast<double>(size) * std_a * stand_dev(b, size, mean_b)), 2.0);
+    }
+
+    /**
+     * @brief The Pearson correlation for two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param size the size of the vector
+     * @param mean_a the mean of the a vector
+     * @param std_a the standard deviation of the a vector
+     * @param mean_b the mean of the b vector
+     * @param std_b the standard deviation of the b vector
+     * @return The coefficient of determination between vector a and vector b
+     */
+    inline double r2(double* a, double* b, int size, double mean_a, double std_a, double mean_b, double std_b)
+    {
+        return std::pow(std::inner_product(a, a + size, b, -1.0 * static_cast<double>(size) * mean_a * mean_b) / (static_cast<double>(size) * std_a * std_b), 2.0);
     }
 
     /**
@@ -183,6 +307,32 @@ namespace util_funcs
      */
     double r2(double* a, double* b, std::vector<int>& sizes);
 
+    /**
+     * @brief Calculate the average Coefficient of Determination between two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param sizes the sizes of the tasks to calculate the correlation on
+     * @param mean_a the mean of the a vector for each task
+     * @param std_a the standard deviation of the a vector for each task
+     * @return The average Coefficient of Determination
+     */
+    double r2(double* a, double* b, std::vector<int>& sizes, std::vector<double>& mean_a, std::vector<double>& std_a);
+
+    /**
+     * @brief Calculate the average Coefficient of Determination between two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param sizes the sizes of the tasks to calculate the correlation on
+     * @param mean_a the mean of the a vector for each task
+     * @param std_a the standard deviation of the a vector for each task
+     * @param mean_b the mean of the b vector for each task
+     * @param std_b the standard deviation of the b vector for each task
+     * @return The average Coefficient of Determination
+     */
+    double r2(double* a, double* b, std::vector<int>& sizes, std::vector<double>& mean_a, std::vector<double>& std_a, std::vector<double>& mean_b, std::vector<double>& std_b);
+
     /**
      * @brief Calculate the average Coefficient of Determination between two vectors
      *
@@ -194,6 +344,34 @@ namespace util_funcs
      */
     double r2(double* a, double* b, int* sz, int n_tasks);
 
+    /**
+     * @brief Calculate the average Coefficient of Determination between two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param sz the start of vector that describes the sizes of the tasks to calculate the correlation on
+     * @param mean_a the mean of the a vector for each task
+     * @param std_a the standard deviation of the a vector for each task
+     * @param n_tasks number of tasks to average over
+     * @return The average Coefficient of Determination
+     */
+    double r2(double* a, double* b, int* sz, double* mean_a, double* std_a, int n_tasks);
+
+    /**
+     * @brief Calculate the average Coefficient of Determination between two vectors
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param sz the start of vector that describes the sizes of the tasks to calculate the correlation on
+     * @param mean_a the mean of the a vector for each task
+     * @param std_a the standard deviation of the a vector for each task
+     * @param mean_b the mean of the b vector for each task
+     * @param std_b the standard deviation of the b vector for each task
+     * @param n_tasks number of tasks to average over
+     * @return The average Coefficient of Determination
+     */
+    double r2(double* a, double* b, int* sz, double* mean_a, double* std_a, double* mean_b, double* std_b, int n_tasks);
+
     /**
      * @brief Calculate the Coefficient of Determination between two vectors (For the log transformed problem)
      *
@@ -205,6 +383,19 @@ namespace util_funcs
      */
     double log_r2(double* a, double* b, double* log_a, int size);
 
+    /**
+     * @brief Calculate the Coefficient of Determination between two vectors (For the log transformed problem)
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param log_a the pointer to the head of the vector used to store the log_transformed a value
+     * @param size the size of the vector
+     * @param mean_b the mean of the b vector for each task
+     * @param std_b the standard deviation of the b vector for each task
+     * @return The Coefficient of Determination
+     */
+    double log_r2(double* a, double* b, double* log_a, int size, double mean_b, double std_b);
+
     /**
      * @brief Calculate the average Coefficient of Determination between two vectors (For the log transformed problem)
      *
@@ -216,6 +407,19 @@ namespace util_funcs
      */
     double log_r2(double* a, double* b, double* log_a, std::vector<int>& sizes);
 
+    /**
+     * @brief Calculate the average Coefficient of Determination between two vectors (For the log transformed problem)
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param log_a the pointer to the head of the vector used to store the log_transformed a value
+     * @param sizes the sizes of the tasks to calculate the correlation on
+     * @param mean_b the mean of the b vector for each task
+     * @param std_b the standard deviation of the b vector for each task
+     * @return The average Coefficient of Determination
+     */
+    double log_r2(double* a, double* b, double* log_a, std::vector<int>& sizes, std::vector<double>& mean_b, std::vector<double>& std_b);
+
     /**
      * @brief Calculate the average Coefficient of Determination between two vectors (For the log transformed problem)
      *
@@ -228,6 +432,20 @@ namespace util_funcs
      */
     double log_r2(double* a, double* b, double* log_a, int* sz, int n_tasks);
 
+    /**
+     * @brief Calculate the average Coefficient of Determination between two vectors (For the log transformed problem)
+     *
+     * @param a the pointer to the head of the first vector
+     * @param b the pointer to the head of the second vector
+     * @param log_a the pointer to the head of the vector used to store the log_transformed a value
+     * @param sz the start of vector that describes the sizes of the tasks to calculate the correlation on
+     * @param mean_b the mean of the b vector for each task
+     * @param std_b the standard deviation of the b vector for each task
+     * @param n_tasks number of tasks to average over
+     * @return The average Coefficient of Determination
+     */
+    double log_r2(double* a, double* b, double* log_a, int* sz, double* mean_b, double* std_b, int n_tasks);
+
     /**
      * @brief Gets the rank variables for a vector
      *
diff --git a/src/utils/project.cpp b/src/utils/project.cpp
index 0e444b88457eb5454e3970bf1fc88b971399b7cb..54e207e1a8799836b9b443daa234169327bcddb3 100644
--- a/src/utils/project.cpp
+++ b/src/utils/project.cpp
@@ -38,14 +38,30 @@ void project_funcs::project_r(double* prop, double* scores, std::vector<node_ptr
         int start = omp_get_thread_num() * (phi.size() / omp_get_num_threads()) + std::min(omp_get_thread_num(), static_cast<int>(phi.size()) % omp_get_num_threads());
         int end = (omp_get_thread_num() + 1) * (phi.size() / omp_get_num_threads()) + std::min(omp_get_thread_num() + 1, static_cast<int>(phi.size()) % omp_get_num_threads());
 
-        std::transform(phi.begin(), phi.end(), scores, [&prop, &sizes](node_ptr feat){return -1 * util_funcs::r(prop, feat->value_ptr(), sizes.data(), sizes.size());});
-
-        for(int pp = 1; pp < n_prop; ++pp)
+        std::vector<double> mean_prop(sizes.size(), 0.0);
+        std::vector<double> std_prop(sizes.size(), 0.0);
+        std::fill_n(scores + start, end - start, 0.0);
+        for(int pp = 0; pp < n_prop; ++pp)
         {
+            int pos = 0;
+            for(int tt = 0; tt < sizes.size(); ++tt)
+            {
+                mean_prop[tt] = util_funcs::mean(prop + pos, sizes[tt]);
+                std_prop[tt] = util_funcs::stand_dev(prop + pos, sizes[tt], mean_prop[tt]);
+                pos += sizes[tt];
+            }
+            std::transform(
+                phi.begin() + start,
+                phi.begin() + end,
+                scores + start,
+                scores + start,
+                [&prop, &sizes, &mean_prop, &std_prop](node_ptr feat, double sc){
+                    return std::min(-1 * std::abs(util_funcs::r(prop, feat->value_ptr(0, true), sizes.data(), mean_prop.data(), std_prop.data(), sizes.size())), sc);
+                }
+            );
             prop += n_samp;
-            std::transform(phi.begin(), phi.end(), scores, scores, [&prop, &sizes](node_ptr feat, double sc){return std::min(-1 * util_funcs::r(prop, feat->value_ptr(), sizes.data(), sizes.size()), sc);});
         }
-        std::transform(scores, scores + phi.size(), scores, [](double score){return std::isnan(score) ? 0.0 : score;});
+        std::transform(scores + start, scores + end, scores + start, [](double score){return std::isnan(score) ? 0.0 : score;});
 
         #pragma omp barrier
     }
@@ -60,12 +76,28 @@ void project_funcs::project_r2(double* prop, double* scores, std::vector<node_pt
         int start = omp_get_thread_num() * (phi.size() / omp_get_num_threads()) + std::min(omp_get_thread_num(), static_cast<int>(phi.size()) % omp_get_num_threads());
         int end = (omp_get_thread_num() + 1) * (phi.size() / omp_get_num_threads()) + std::min(omp_get_thread_num() + 1, static_cast<int>(phi.size()) % omp_get_num_threads());
 
-        std::transform(phi.begin() + start, phi.begin() + end, scores + start, [&prop, &sizes](node_ptr feat){return -1 * util_funcs::r2(prop, feat->value_ptr(), sizes.data(), sizes.size());});
-
-        for(int pp = 1; pp < n_prop; ++pp)
+        std::vector<double> mean_prop(sizes.size(), 0.0);
+        std::vector<double> std_prop(sizes.size(), 0.0);
+        std::fill_n(scores + start, end - start, 0.0);
+        for(int pp = 0; pp < n_prop; ++pp)
         {
+            int pos = 0;
+            for(int tt = 0; tt < sizes.size(); ++tt)
+            {
+                mean_prop[tt] = util_funcs::mean(prop + pos, sizes[tt]);
+                std_prop[tt] = util_funcs::stand_dev(prop + pos, sizes[tt], mean_prop[tt]);
+                pos += sizes[tt];
+            }
+            std::transform(
+                phi.begin() + start,
+                phi.begin() + end,
+                scores + start,
+                scores + start,
+                [&prop, &sizes, &mean_prop, &std_prop](node_ptr feat, double sc){
+                    return std::min(-1 * util_funcs::r2(prop, feat->value_ptr(0, true), sizes.data(), mean_prop.data(), std_prop.data(), sizes.size()), sc);
+                }
+            );
             prop += n_samp;
-            std::transform(phi.begin() + start, phi.begin() + end, scores + start, scores + start, [&prop, &sizes](node_ptr feat, double sc){return std::min(-1 * util_funcs::r2(prop, feat->value_ptr(), sizes.data(), sizes.size()), sc);});
         }
         std::transform(scores + start, scores + end, scores + start, [](double score){return std::isnan(score) ? 0.0 : score;});
 
@@ -82,15 +114,33 @@ void project_funcs::project_log_r2(double* prop, double* scores, std::vector<nod
         int start = omp_get_thread_num() * (phi.size() / omp_get_num_threads()) + std::min(omp_get_thread_num(), static_cast<int>(phi.size()) % omp_get_num_threads());
         int end = (omp_get_thread_num() + 1) * (phi.size() / omp_get_num_threads()) + std::min(omp_get_thread_num() + 1, static_cast<int>(phi.size()) % omp_get_num_threads());
 
+        std::vector<double> mean_prop(sizes.size(), 0.0);
+        std::vector<double> std_prop(sizes.size(), 0.0);
         std::vector<double> log_feat(n_samp, 0.0);
-        std::transform(phi.begin() + start, phi.begin() + end, scores + start, [&prop, &sizes, &log_feat](node_ptr feat){return -1 * util_funcs::log_r2(feat->value_ptr(), prop, log_feat.data(), sizes.data(), sizes.size());});
 
-        for(int pp = 1; pp < n_prop; ++pp)
+        std::fill_n(scores + start, end - start, 0.0);
+        for(int pp = 0; pp < n_prop; ++pp)
         {
+            int pos = 0;
+            for(int tt = 0; tt < sizes.size(); ++tt)
+            {
+                mean_prop[tt] = util_funcs::mean(prop + pos, sizes[tt]);
+                std_prop[tt] = util_funcs::stand_dev(prop + pos, sizes[tt], mean_prop[tt]);
+                pos += sizes[tt];
+            }
+            std::transform(
+                phi.begin() + start,
+                phi.begin() + end,
+                scores + start,
+                scores + start,
+                [&prop, &sizes, &mean_prop, &std_prop, &log_feat](node_ptr feat, double sc){
+                    return std::min(-1 * util_funcs::log_r2(feat->value_ptr(0, true), prop, log_feat.data(), sizes.data(), mean_prop.data(), std_prop.data(), sizes.size()), sc);
+                }
+            );
             prop += n_samp;
-            std::transform(phi.begin() + start, phi.begin() + end, scores + start, scores + start, [&prop, &sizes, &log_feat](node_ptr feat, double sc){return std::min(-1 * util_funcs::log_r2(feat->value_ptr(), prop, log_feat.data(), sizes.data(), sizes.size()), sc);});
         }
         std::transform(scores + start, scores + end, scores + start, [](double score){return std::isnan(score) ? 0.0 : score;});
+
         #pragma omp barrier
     }
 }
@@ -107,7 +157,7 @@ void project_funcs::project_classify(double* prop, double* scores, std::vector<n
         for(int pp = 0; pp < n_prop; ++pp)
         {
             ConvexHull1D convex_hull(sizes, prop);
-            std::transform(phi.begin() + start, phi.begin() + end, scores + start, scores + start, [&convex_hull](node_ptr feat, double score){return std::min(convex_hull.overlap_1d(feat->value_ptr()), score);});
+            std::transform(phi.begin() + start, phi.begin() + end, scores + start, scores + start, [&convex_hull](node_ptr feat, double score){return std::min(convex_hull.overlap_1d(feat->value_ptr(0, true)), score);});
             prop += n_samp;
             #pragma omp barrier
         }
@@ -119,12 +169,30 @@ void project_funcs::project_r_no_omp(double* prop, double* scores, std::vector<n
 {
     int n_samp = std::accumulate(sizes.begin(), sizes.end(), 0);
 
-    std::transform(phi.begin(), phi.end(), scores, [&prop, &sizes](node_ptr feat){return -1 * util_funcs::r(prop, feat->value_ptr(), sizes.data(), sizes.size());});
+    std::vector<double> mean_prop(sizes.size(), 0.0);
+    std::vector<double> std_prop(sizes.size(), 0.0);
 
-    for(int pp = 1; pp < n_prop; ++pp)
+    std::fill_n(scores, phi.size(), 0.0);
+
+    for(int pp = 0; pp < n_prop; ++pp)
     {
+        int pos = 0;
+        for(int tt = 0; tt < sizes.size(); ++tt)
+        {
+            mean_prop[tt] = util_funcs::mean(prop + pos, sizes[tt]);
+            std_prop[tt] = util_funcs::stand_dev(prop + pos, sizes[tt], mean_prop[tt]);
+            pos += sizes[tt];
+        }
+        std::transform(
+            phi.begin(),
+            phi.end(),
+            scores,
+            scores,
+            [&prop, &sizes, &mean_prop, &std_prop](node_ptr feat, double sc){
+                return std::min(-1 * util_funcs::r(prop, feat->value_ptr(0, true), sizes.data(), mean_prop.data(), std_prop.data(), sizes.size()), sc);
+            }
+        );
         prop += n_samp;
-        std::transform(phi.begin(), phi.end(), scores, scores, [&prop, &sizes](node_ptr feat, double sc){return std::min(-1 * util_funcs::r(prop, feat->value_ptr(), sizes.data(), sizes.size()), sc);});
     }
     std::transform(scores, scores + phi.size(), scores, [](double score){return std::isnan(score) ? 0.0 : score;});
 }
@@ -133,12 +201,30 @@ void project_funcs::project_r2_no_omp(double* prop, double* scores, std::vector<
 {
     int n_samp = std::accumulate(sizes.begin(), sizes.end(), 0);
 
-    std::transform(phi.begin(), phi.end(), scores, [&prop, &sizes](node_ptr feat){return -1 * util_funcs::r2(prop, feat->value_ptr(), sizes.data(), sizes.size());});
+    std::vector<double> mean_prop(sizes.size(), 0.0);
+    std::vector<double> std_prop(sizes.size(), 0.0);
+
+    std::fill_n(scores, phi.size(), 0.0);
 
-    for(int pp = 1; pp < n_prop; ++pp)
+    for(int pp = 0; pp < n_prop; ++pp)
     {
+        int pos = 0;
+        for(int tt = 0; tt < sizes.size(); ++tt)
+        {
+            mean_prop[tt] = util_funcs::mean(prop + pos, sizes[tt]);
+            std_prop[tt] = util_funcs::stand_dev(prop + pos, sizes[tt], mean_prop[tt]);
+            pos += sizes[tt];
+        }
+        std::transform(
+            phi.begin(),
+            phi.end(),
+            scores,
+            scores,
+            [&prop, &sizes, &mean_prop, &std_prop](node_ptr feat, double sc){
+                return std::min(-1 * util_funcs::r2(prop, feat->value_ptr(0, true), sizes.data(), mean_prop.data(), std_prop.data(), sizes.size()), sc);
+            }
+        );
         prop += n_samp;
-        std::transform(phi.begin(), phi.end(), scores, scores, [&prop, &sizes](node_ptr feat, double sc){return std::min(-1 * util_funcs::r2(prop, feat->value_ptr(), sizes.data(), sizes.size()), sc);});
     }
     std::transform(scores, scores + phi.size(), scores, [](double score){return std::isnan(score) ? 0.0 : score;});
 }
@@ -148,13 +234,30 @@ void project_funcs::project_log_r2_no_omp(double* prop, double* scores, std::vec
     int n_samp = std::accumulate(sizes.begin(), sizes.end(), 0);
 
     std::vector<double> log_feat(n_samp, 0.0);
+    std::vector<double> mean_prop(sizes.size(), 0.0);
+    std::vector<double> std_prop(sizes.size(), 0.0);
 
-    std::transform(phi.begin(), phi.end(), scores, [&prop, &sizes, &log_feat](node_ptr feat){return -1 * util_funcs::log_r2(feat->value_ptr(), prop, log_feat.data(), sizes.data(), sizes.size());});
+    std::fill_n(scores, phi.size(), 0.0);
 
-    for(int pp = 1; pp < n_prop; ++pp)
+    for(int pp = 0; pp < n_prop; ++pp)
     {
+        int pos = 0;
+        for(int tt = 0; tt < sizes.size(); ++tt)
+        {
+            mean_prop[tt] = util_funcs::mean(prop + pos, sizes[tt]);
+            std_prop[tt] = util_funcs::stand_dev(prop + pos, sizes[tt], mean_prop[tt]);
+            pos += sizes[tt];
+        }
+        std::transform(
+            phi.begin(),
+            phi.end(),
+            scores,
+            scores,
+            [&prop, &sizes, &mean_prop, &std_prop, &log_feat](node_ptr feat, double sc){
+                return std::min(-1 * util_funcs::log_r2(feat->value_ptr(0, true), prop, log_feat.data(), sizes.data(), mean_prop.data(), std_prop.data(), sizes.size()), sc);
+            }
+        );
         prop += n_samp;
-        std::transform(phi.begin(), phi.end(), scores, scores, [&prop, &sizes, &log_feat](node_ptr feat, double sc){return std::min(-1 * util_funcs::log_r2(feat->value_ptr(), prop, log_feat.data(), sizes.data(), sizes.size()), sc);});
     }
     std::transform(scores, scores + phi.size(), scores, [](double score){return std::isnan(score) ? 0.0 : score;});
 }
@@ -167,7 +270,7 @@ void project_funcs::project_classify_no_omp(double* prop, double* scores, std::v
     for(int pp = 0; pp < n_prop; ++pp)
     {
         ConvexHull1D convex_hull(sizes, prop);
-        std::transform(phi.begin(), phi.end(), scores, scores, [&convex_hull](node_ptr feat, double score){return std::min(convex_hull.overlap_1d(feat->value_ptr()), score);});
+        std::transform(phi.begin(), phi.end(), scores, scores, [&convex_hull](node_ptr feat, double score){return std::min(convex_hull.overlap_1d(feat->value_ptr(0, true)), score);});
         prop += n_samp;
     }
     std::transform(scores, scores + phi.size(), scores, [](double score){return std::isnan(score) ? std::numeric_limits<double>::infinity() : score;});
diff --git a/tests/googletest/feature_creation/value_storage/test_value_storage.cc b/tests/googletest/feature_creation/value_storage/test_value_storage.cc
index 13f7af5541c0809562857746d0ef8a5524a2a1b8..d866e473df58584bf431d62ef59bbad4451ab38b 100644
--- a/tests/googletest/feature_creation/value_storage/test_value_storage.cc
+++ b/tests/googletest/feature_creation/value_storage/test_value_storage.cc
@@ -12,14 +12,14 @@ namespace {
         EXPECT_EQ(node_value_arrs::N_SAMPLES_TEST, 2);
         EXPECT_EQ(node_value_arrs::N_RUNGS_STORED, 0);
         EXPECT_EQ(node_value_arrs::N_STORE_FEATURES, 1);
-        EXPECT_EQ(node_value_arrs::N_OP_SLOTS, 3);
+        EXPECT_EQ(node_value_arrs::N_OP_SLOTS, 6);
         EXPECT_EQ(node_value_arrs::MAX_RUNG, 2);
         EXPECT_EQ(node_value_arrs::VALUES_ARR.size(), 5);
         EXPECT_EQ(node_value_arrs::TEST_VALUES_ARR.size(), 2);
-        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_ARR.size(), node_value_arrs::MAX_N_THREADS * (3 * 1 + 1) * 5);
-        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_REG.size(), node_value_arrs::MAX_N_THREADS * (3 * 1 + 1));
-        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_TEST_ARR.size(), node_value_arrs::MAX_N_THREADS * (3 * 1 + 1) * 2);
-        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_TEST_REG.size(), node_value_arrs::MAX_N_THREADS * (3 * 1 + 1));
+        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_ARR.size(), node_value_arrs::MAX_N_THREADS * (6 * 1 + 1) * 5);
+        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_REG.size(), node_value_arrs::MAX_N_THREADS * (6 * 1 + 1));
+        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_TEST_ARR.size(), node_value_arrs::MAX_N_THREADS * (6 * 1 + 1) * 2);
+        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_TEST_REG.size(), node_value_arrs::MAX_N_THREADS * (6 * 1 + 1));
 
         node_value_arrs::resize_values_arr(1, 2, false);
         EXPECT_EQ(node_value_arrs::N_SAMPLES, 5);
@@ -39,13 +39,13 @@ namespace {
         EXPECT_EQ(node_value_arrs::N_SAMPLES_TEST, 2);
         EXPECT_EQ(node_value_arrs::N_RUNGS_STORED, 1);
         EXPECT_EQ(node_value_arrs::N_STORE_FEATURES, 2);
-        EXPECT_EQ(node_value_arrs::N_OP_SLOTS, 1);
+        EXPECT_EQ(node_value_arrs::N_OP_SLOTS, 2);
         EXPECT_EQ(node_value_arrs::VALUES_ARR.size(), 10);
         EXPECT_EQ(node_value_arrs::TEST_VALUES_ARR.size(), 4);
-        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_ARR.size(), node_value_arrs::MAX_N_THREADS * (1 * 2 + 1) * 5);
-        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_REG.size(), node_value_arrs::MAX_N_THREADS * (1 * 2 + 1));
-        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_TEST_ARR.size(), node_value_arrs::MAX_N_THREADS * (1 * 2 + 1) * 2);
-        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_TEST_REG.size(), node_value_arrs::MAX_N_THREADS * (1 * 2 + 1));
+        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_ARR.size(), node_value_arrs::MAX_N_THREADS * (2 * 2 + 1) * 5);
+        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_REG.size(), node_value_arrs::MAX_N_THREADS * (2 * 2 + 1));
+        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_TEST_ARR.size(), node_value_arrs::MAX_N_THREADS * (2 * 2 + 1) * 2);
+        EXPECT_EQ(node_value_arrs::TEMP_STORAGE_TEST_REG.size(), node_value_arrs::MAX_N_THREADS * (2 * 2 + 1));
 
         node_value_arrs::initialize_d_matrix_arr();
         EXPECT_EQ(node_value_arrs::N_SELECTED, 0);
@@ -66,11 +66,11 @@ namespace {
         EXPECT_EQ(node_value_arrs::TEST_VALUES_ARR[3], 1.0);
 
         node_value_arrs::get_value_ptr(10, 141, 2, 0)[1] = 1.0;
-        EXPECT_EQ(node_value_arrs::temp_storage_reg(10, 2, 0), 141);
+        EXPECT_EQ(node_value_arrs::temp_storage_reg(10, 2, 0, false), 141);
         EXPECT_EQ(node_value_arrs::access_temp_storage((10 % 2) + omp_get_thread_num() * (2 * 1 + 1))[1], 1.0);
 
         node_value_arrs::get_test_value_ptr(10, 141, 2, 0)[1] = 1.0;
-        EXPECT_EQ(node_value_arrs::temp_storage_test_reg(10, 2, 0), 141);
+        EXPECT_EQ(node_value_arrs::temp_storage_test_reg(10, 2, 0, false), 141);
         EXPECT_EQ(node_value_arrs::access_temp_storage_test((10 % 2) + omp_get_thread_num() * (2 * 1 + 1))[1], 1.0);
 
         node_value_arrs::get_d_matrix_ptr(1)[0] = 1.0;
@@ -78,10 +78,10 @@ namespace {
 
         #pragma omp parallel
         {
-            std::fill_n(node_value_arrs::TEMP_STORAGE_REG.data() + 2 * omp_get_thread_num(), 2, omp_get_thread_num() + 1);
-            EXPECT_EQ(node_value_arrs::TEMP_STORAGE_REG[2 * omp_get_thread_num()], omp_get_thread_num() + 1);
+            std::fill_n(node_value_arrs::TEMP_STORAGE_REG.data() + 5 * omp_get_thread_num(), 5, omp_get_thread_num() + 1);
+            EXPECT_EQ(node_value_arrs::TEMP_STORAGE_REG[5 * omp_get_thread_num()], omp_get_thread_num() + 1);
             node_value_arrs::clear_temp_reg_thread();
-            EXPECT_EQ(node_value_arrs::TEMP_STORAGE_REG[2 * omp_get_thread_num()], -1);
+            EXPECT_EQ(node_value_arrs::TEMP_STORAGE_REG[5 * omp_get_thread_num()], -1);
         }
 
         std::fill_n(node_value_arrs::TEMP_STORAGE_REG.data(), node_value_arrs::TEMP_STORAGE_REG.size(), 2.0);