4 #include <parcer/Eigen.h>
7 namespace pt = boost::property_tree;
11 LocalRegression::LocalRegression(std::shared_ptr<ModPiece>
function, pt::ptree& pt) :
ModPiece(function->inputSizes, function->outputSizes), kn(pt.
get<unsigned int>(
"NumNeighbors")) {
16 LocalRegression::LocalRegression(std::shared_ptr<muq::Modeling::ModPiece>
function, boost::property_tree::ptree& pt, std::shared_ptr<parcer::Communicator> comm) :
ModPiece(function->inputSizes, function->outputSizes), kn(pt.
get<unsigned int>(
"NumNeighbors")), comm(comm) {
37 cache = std::make_shared<FlannCache>(
function);
40 pt.put<std::string>(
"PolynomialBasis", pt.get<std::string>(
"PolynomialBasis",
"Legendre"));
41 pt.put<
unsigned int>(
"Order", pt.get<
unsigned int>(
"Order", 2));
42 pt.put<
double>(
"MaxPoisednessRadius", pt.get<
double>(
"MaxPoisednessRadius", 1.0));
43 pt.put<
unsigned int>(
"InputSize",
function->inputSizes(0));
44 reg = std::make_shared<Regression>(pt);
49 cache = std::make_shared<FlannCache>(
cache->Function());
54 std::vector<Eigen::VectorXd> neighbors;
55 std::vector<Eigen::VectorXd> result;
56 cache->NearestNeighbors(input,
kn, neighbors, result);
59 reg->Fit(neighbors, result, input);
72 outputs[0] = (Eigen::VectorXd)boost::any_cast<Eigen::MatrixXd const&>(
reg->Evaluate(inputs[0].get()) [0]).col(0);
82 return cache->at(index);
87 return cache->InCache(point)>=0;
92 const Eigen::VectorXd result =
cache->Add(input);
96 for(
unsigned int i=0; i<
comm->GetSize(); ++i ) {
97 if( i==
comm->GetRank() ) {
continue; }
99 comm->Send(std::pair<Eigen::VectorXd, Eigen::VectorXd>(input, result), i,
tagSingle);
109 for(
auto it : inputs ) {
Add(it); }
114 std::vector<Eigen::VectorXd> neighbors;
115 cache->NearestNeighbors(input,
kn, neighbors);
121 assert(neighbors.size()>=
kn);
123 std::pair<Eigen::VectorXd, double> lambda =
reg->PoisednessConstant(neighbors, input,
kn);
126 unsigned int index = 0;
127 for(
unsigned int i=0; i<
kn; ++i ) {
128 const double d = (input-neighbors[i]).norm();
129 if( d<dist ) { dist = d; index = i; }
132 return std::tuple<Eigen::VectorXd, double, unsigned int>(lambda.first, lambda.second, index);
137 std::vector<Eigen::VectorXd> neighbors;
138 cache->NearestNeighbors(input,
kn, neighbors);
146 assert(neighbors.size()==
kn);
153 for(
auto n : neighbors) { delta = std::max(delta, (
n-input).norm()); }
160 return std::pair<double, double>(((
double)
reg->order+1.0)*std::log(delta), delta);
165 cache->NearestNeighbors(input,
kn, neighbors);
170 cache->NearestNeighbors(input,
kn, neighbors, result);
179 reg->Fit(neighbors, result, input);
182 return (Eigen::VectorXd)boost::any_cast<Eigen::MatrixXd const&>(
reg->Evaluate(input) [0]).col(0);
187 if( !
comm ) {
return; }
190 for(
unsigned int i=0; i<
comm->GetSize(); ++i ) {
193 if( i==
comm->GetRank() ) {
continue; }
200 parcer::RecvRequest recvReq;
207 const std::pair<Eigen::VectorXd, Eigen::VectorXd> point = recvReq.GetObject<std::pair<Eigen::VectorXd, Eigen::VectorXd> >();
214 cache->Add(point.first, point.second);
238 return cache->Centroid();
std::pair< double, double > ErrorIndicator(Eigen::VectorXd const &input) const
Get the error indicator.
Eigen::VectorXd CacheCentroid() const
Get the centroid of the cache.
virtual void EvaluateImpl(muq::Modeling::ref_vector< Eigen::VectorXd > const &inputs) override
void SetUp(std::shared_ptr< muq::Modeling::ModPiece > function, boost::property_tree::ptree &pt)
Set up the regressor.
unsigned int Order() const
Polynomial order.
Eigen::VectorXd EvaluateRegressor(Eigen::VectorXd const &input, std::vector< Eigen::VectorXd > const &neighbors, std::vector< Eigen::VectorXd > const &result) const
void ClearCache()
Clear the cache.
Eigen::VectorXd CachePoint(unsigned int const index) const
A point in the cache.
unsigned int CacheSize() const
Get the total size of the cache.
std::shared_ptr< muq::Modeling::FlannCache > cache
A cache containing previous model evaluations.
std::tuple< Eigen::VectorXd, double, unsigned int > PoisednessConstant(Eigen::VectorXd const &input) const
Get the poisedness constant.
std::shared_ptr< parcer::Communicator > comm
void Add(std::vector< Eigen::VectorXd > const &inputs) const
Add some points to the cache.
std::shared_ptr< Regression > reg
A regressor.
const unsigned int kn
The number of nearest neighbors to use by the regressor.
bool InCache(Eigen::VectorXd const &point) const
Is a point in the cache?
void FitRegression(Eigen::VectorXd const &input) const
Fit the regression to the nearest neighbors.
void NearestNeighbors(Eigen::VectorXd const &input, std::vector< Eigen::VectorXd > &neighbors) const
Get the number of nearest neighbors.
Provides an abstract interface for defining vector-valued model components.
const Eigen::VectorXi inputSizes
std::vector< Eigen::VectorXd > outputs
const Eigen::VectorXi outputSizes
std::vector< std::reference_wrapper< const T > > ref_vector
A vector of references to something ...
auto get(const nlohmann::detail::iteration_proxy_value< IteratorType > &i) -> decltype(i.key())