.. _program_listing_file_include_gwmodelpp_GWRMultiscale.h:

Program Listing for File GWRMultiscale.h
========================================

|exhale_lsh| :ref:`Return to documentation for file <file_include_gwmodelpp_GWRMultiscale.h>` (``include/gwmodelpp/GWRMultiscale.h``)

.. |exhale_lsh| unicode:: U+021B0 .. UPWARDS ARROW WITH TIP LEFTWARDS

.. code-block:: cpp

   #ifndef GWRMULTISCALE_H
   #define GWRMULTISCALE_H
   
   #include <utility>
   #include <string>
   #include <initializer_list>
   #include <optional>
   #include "SpatialMultiscaleAlgorithm.h"
   #include "spatialweight/SpatialWeight.h"
   #include "IRegressionAnalysis.h"
   #include "IBandwidthSelectable.h"
   #include "BandwidthSelector.h"
   #include "IParallelizable.h"
   
   
   namespace gwm
   {
   
   #define GWM_LOG_TAG_MGWR_INITIAL_BW "#initial-bandwidth "
   #define GWM_LOG_TAG_MGWR_BACKFITTING "#back-fitting "
   
   #define GWM_LOG_MGWR_BACKFITTING(MESSAGE) { GWM_LOG_INFO((std::string(GWM_LOG_TAG_MGWR_BACKFITTING) + (MESSAGE))); }
   
   class GWRMultiscale : public SpatialMultiscaleAlgorithm,
       public IBandwidthSelectable,
       public IRegressionAnalysis,
       public IParallelizable,
       public IParallelOpenmpEnabled,
       public IParallelCudaEnabled,
       public IParallelMpiEnabled
   {
   public:
   
       enum BandwidthInitilizeType
       {
           Null,
           Initial,
           Specified
       };
   
       static std::unordered_map<BandwidthInitilizeType,std::string> BandwidthInitilizeTypeNameMapper; 
   
       enum BandwidthSelectionCriterionType
       {
           CV,
           AIC
       };
   
       static std::unordered_map<BandwidthSelectionCriterionType,std::string> BandwidthSelectionCriterionTypeNameMapper;   
   
       enum BackFittingCriterionType
       {
           CVR,
           dCVR
       };
   
       static std::unordered_map<BackFittingCriterionType,std::string> BackFittingCriterionTypeNameMapper; 
   
       typedef double (GWRMultiscale::*BandwidthSizeCriterionFunction)(BandwidthWeight*);  
           
       typedef arma::vec (GWRMultiscale::*FitVarFunction)(const size_t);  
           
       typedef arma::vec (GWRMultiscale::*FitVarCoreFunction)(const arma::vec&, const arma::vec&, const SpatialWeight&, arma::mat&);  
           
       typedef arma::vec (GWRMultiscale::*FitVarCoreCVFunction)(const arma::vec&, const arma::vec&, const SpatialWeight&);  
           
       typedef arma::vec (GWRMultiscale::*FitVarCoreSHatFunction)(const arma::vec&, const arma::vec&, const SpatialWeight&, arma::vec&);  
       
   private:
   
       static arma::vec Fitted(const arma::mat& x, const arma::mat& betas)
       {
           return sum(betas % x, 1);
       }
   
       static double RSS(const arma::mat& x, const arma::mat& y, const arma::mat& betas)
       {
           arma::vec r = y - Fitted(x, betas);
           return sum(r % r);
       }
   
       static double AICc(const arma::mat& x, const arma::mat& y, const arma::mat& betas, const arma::vec& shat)
       {
           double ss = RSS(x, y, betas), n = (double)x.n_rows;
           return n * log(ss / n) + n * log(2 * arma::datum::pi) + n * ((n + shat(0)) / (n - 2 - shat(0)));
       }
   
       bool isStoreS()
       {
           return mHasHatMatrix && (mCoords.n_rows < 8192);
       }
   
       static RegressionDiagnostic CalcDiagnostic(const arma::mat &x, const arma::vec &y, const arma::vec &shat, double RSS);
   
       void setInitSpatialWeight(const SpatialWeight &spatialWeight)
       {
           mInitSpatialWeight = spatialWeight;
       }
   
   public:
   
       GWRMultiscale() {}
   
       GWRMultiscale(const arma::mat& x, const arma::vec& y, const arma::mat& coords, const std::vector<SpatialWeight>& spatialWeights)
           : SpatialMultiscaleAlgorithm(coords, spatialWeights)
       {
           mX = x;
           mY = y;
           if (spatialWeights.size() > 0)
               setInitSpatialWeight(spatialWeights[0]);
       }
   
       virtual ~GWRMultiscale() {}
   
   public:
   
       void setGoldenUpperBounds(double value) { mGoldenUpperBounds = value; }
   
       void setGoldenLowerBounds(double value) { mGoldenLowerBounds = value; }
   
       const std::vector<BandwidthInitilizeType>& bandwidthInitilize() const { return GWRMultiscale::mBandwidthInitilize; }
       
       void setBandwidthInitilize(const std::vector<BandwidthInitilizeType> &bandwidthInitilize);
   
       const std::vector<BandwidthSelectionCriterionType>& bandwidthSelectionApproach() const { return GWRMultiscale::mBandwidthSelectionApproach; }
       
       void setBandwidthSelectionApproach(const std::vector<BandwidthSelectionCriterionType> &bandwidthSelectionApproach);
   
       const std::vector<bool>& preditorCentered() const { return mPreditorCentered; }
       
       void setPreditorCentered(const std::vector<bool> &preditorCentered) { mPreditorCentered = preditorCentered; }
   
       const std::vector<double>& bandwidthSelectThreshold() const { return mBandwidthSelectThreshold; }
       
       void setBandwidthSelectThreshold(const std::vector<double> &bandwidthSelectThreshold) { mBandwidthSelectThreshold = bandwidthSelectThreshold; }
   
       bool hasHatMatrix() const { return mHasHatMatrix; }
       
       void setHasHatMatrix(bool hasHatMatrix) { mHasHatMatrix = hasHatMatrix; }
   
       size_t bandwidthSelectRetryTimes() const { return (size_t)mBandwidthSelectRetryTimes; }
       
       void setBandwidthSelectRetryTimes(size_t bandwidthSelectRetryTimes) { mBandwidthSelectRetryTimes = (arma::uword)bandwidthSelectRetryTimes; }
   
       size_t maxIteration() const { return mMaxIteration; }
       
       void setMaxIteration(size_t maxIteration) { mMaxIteration = maxIteration; }
   
       BackFittingCriterionType criterionType() const { return mCriterionType; }
       
       void setCriterionType(const BackFittingCriterionType &criterionType) { mCriterionType = criterionType; }
   
       double criterionThreshold() const { return mCriterionThreshold; }
       
       void setCriterionThreshold(double criterionThreshold) { mCriterionThreshold = criterionThreshold; }
   
       int adaptiveLower() const { return mAdaptiveLower; }
       
       void setAdaptiveLower(int adaptiveLower) { mAdaptiveLower = adaptiveLower; }
   
       const arma::mat& betas() const { return mBetas; }
       
       BandwidthSizeCriterionFunction bandwidthSizeCriterionVar(BandwidthSelectionCriterionType type);
   
   public:     // SpatialAlgorithm interface
       bool isValid() override;
   
   
   public:     // SpatialMultiscaleAlgorithm interface
       virtual void setSpatialWeights(const std::vector<SpatialWeight> &spatialWeights) override;
   
   
   public:     // IBandwidthSizeSelectable interface
       Status getCriterion(BandwidthWeight* weight, double& criterion) override
       {
           criterion = (this->*mBandwidthSizeCriterion)(weight);
           return mStatus;
       }
   
   
   public:     // IRegressionAnalysis interface
       virtual const arma::vec& dependentVariable() const override { return mY; }
       virtual void setDependentVariable(const arma::vec& y) override { mY = y; }
   
       virtual const arma::mat& independentVariables() const override { return mX; }
       virtual void setIndependentVariables(const arma::mat& x) override { mX = x; }
   
       virtual bool hasIntercept() const override { return mHasIntercept; }
       virtual void setHasIntercept(const bool has) override { mHasIntercept = has; }
   
       virtual RegressionDiagnostic diagnostic() const override { return mDiagnostic; }
   
       arma::mat predict(const arma::mat& locations) override { return arma::mat(locations.n_rows, mX.n_cols, arma::fill::zeros); }
   
       arma::mat fit() override;
   
   public:     // IParallelalbe interface
       int parallelAbility() const override 
       {
           return ParallelType::SerialOnly
   #ifdef ENABLE_OPENMP
               | ParallelType::OpenMP
   #endif
   #ifdef ENABLE_CUDA
               | ParallelType::CUDA
   #endif
               ;
       }
   
       ParallelType parallelType() const override { return mParallelType; }
   
       void setParallelType(const ParallelType &type) override;
   
   
   public:     // IOpenmpParallelable interface
       void setOmpThreadNum(const int threadNum) override { mOmpThreadNum = threadNum; }
   
   public:     // IParallelCudaEnabled interface
       void setGPUId(const int gpuId) override { mGpuId = gpuId; }
       void setGroupSize(const std::size_t size) override { mGroupLength = size; }
       int workerId() override { return mWorkerId; }
       void setWorkerId(int id) override { mWorkerId = id; };
       void setWorkerNum(int size) override { mWorkerNum = size; };
   
   protected:
   
       BandwidthWeight* bandwidth(size_t i)
       {
           return mSpatialWeights[i].weight<BandwidthWeight>();
       }
   
       arma::mat fitInitial();
   
       arma::mat backfitting(const arma::mat& betas0);
   
       arma::vec fitVarBase(const size_t var);
       
       double bandwidthSizeCriterionVarCVBase(BandwidthWeight* bandwidthWeight);
       
       double bandwidthSizeCriterionVarAICBase(BandwidthWeight* bandwidthWeight);
   
       arma::vec fitVarCoreSerial(const arma::vec& x, const arma::vec& y, const SpatialWeight& sw, arma::mat& S);
   
       arma::vec fitVarCoreCVSerial(const arma::vec& x, const arma::vec& y, const SpatialWeight& sw);
   
       arma::vec fitVarCoreSHatSerial(const arma::vec& x, const arma::vec& y, const SpatialWeight& sw, arma::vec& shat);
   
   #ifdef ENABLE_OPENMP
       arma::vec fitVarCoreOmp(const arma::vec& x, const arma::vec& y, const SpatialWeight& sw, arma::mat& S);
   
       arma::vec fitVarCoreCVOmp(const arma::vec& x, const arma::vec& y, const SpatialWeight& sw);
   
       arma::vec fitVarCoreSHatOmp(const arma::vec& x, const arma::vec& y, const SpatialWeight& sw, arma::vec& shat);
   #endif
   
   #ifdef ENABLE_CUDA
       
       arma::vec fitVarCoreCuda(const arma::vec& x, const arma::vec& y, const SpatialWeight& sw, arma::mat& S);
       
       arma::vec fitVarCoreCVCuda(const arma::vec& x, const arma::vec& y, const SpatialWeight& sw);
       
       arma::vec fitVarCoreSHatCuda(const arma::vec& x, const arma::vec& y, const SpatialWeight& sw, arma::vec& shat);
   
   #endif // ENABLE_CUDA
   
   #ifdef ENABLE_MPI
       arma::vec fitVarMpi(const size_t var);
       double bandwidthSizeCriterionVarCVMpi(BandwidthWeight* bandwidthWeight);
       double bandwidthSizeCriterionVarAICMpi(BandwidthWeight* bandwidthWeight);
   #endif // ENABLE_MPI
   
   private:
       FitVarFunction mFitVar = &GWRMultiscale::fitVarBase;  
       FitVarCoreFunction mFitVarCore = &GWRMultiscale::fitVarCoreSerial;  
       FitVarCoreCVFunction mFitVarCoreCV = &GWRMultiscale::fitVarCoreCVSerial;  
       FitVarCoreSHatFunction mFitVarCoreSHat = &GWRMultiscale::fitVarCoreSHatSerial;  
   
       SpatialWeight mInitSpatialWeight;   
       BandwidthSizeCriterionFunction mBandwidthSizeCriterion = &GWRMultiscale::bandwidthSizeCriterionVarCVBase; 
       size_t mBandwidthSelectionCurrentIndex = 0; 
       double mBandwidthLastCriterion = DBL_MAX;   
       std::optional<double> mGoldenUpperBounds;
       std::optional<double> mGoldenLowerBounds;
       std::vector<double> mMaxDistances;
       std::vector<double> mMinDistances;
   
       std::vector<BandwidthInitilizeType> mBandwidthInitilize;    
       std::vector<BandwidthSelectionCriterionType> mBandwidthSelectionApproach;   
       std::vector<bool> mPreditorCentered;    
       std::vector<double> mBandwidthSelectThreshold;  
       arma::uword mBandwidthSelectRetryTimes = 5; 
       size_t mMaxIteration = 500; 
       BackFittingCriterionType mCriterionType = BackFittingCriterionType::dCVR;   
       double mCriterionThreshold = 1e-6;  
       int mAdaptiveLower = 10;    
   
       bool mHasHatMatrix = true;  
   
       arma::mat mX;   
       arma::vec mY;   
       arma::mat mBetas;   
       arma::mat mBetasSE; 
       arma::mat mBetasTV; 
       bool mHasIntercept = true;  
   
       arma::mat mS0;  
       arma::cube mSArray; 
       arma::cube mC;  
       arma::mat mX0;  
       arma::vec mY0;  
       arma::vec mXi;  
       arma::vec mYi;  
   
       double mRSS0;   
   
       RegressionDiagnostic mDiagnostic;   
   
       ParallelType mParallelType = ParallelType::SerialOnly;  
       int mOmpThreadNum = 8;  
       size_t mGroupLength = 64;   
       int mGpuId = 0; 
       int mWorkerId = 0;
       int mWorkerNum = 1;
   #ifdef ENABLE_MPI
       arma::uword mWorkRangeSize = 0;
   #endif // ENABLE_MPI
       std::pair<arma::uword, arma::uword> mWorkRange = std::make_pair(arma::uword(0), arma::uword(0));
   
   public:
       static int treeChildCount;  
   };
   
   }
   
   #endif // GWRMULTISCALE_H