00001
00002
00003
00004
00005 #ifndef MCLIST_H
00006 #define MCLIST_H
00007
00008 #include <cmath>
00009 #include <vector>
00010 #include <algorithm>
00011 #include <iostream>
00012 #include "errors.h"
00013 #include "special.h"
00014 #include "data.h"
00015 #include "rng.h"
00016
00021 class PsiMClist
00022 {
00023 private:
00024 std::vector< std::vector<double> > mcestimates;
00025 std::vector<double> deviances;
00026 public:
00027 PsiMClist (
00028 int N,
00029 int nprm
00030 ) : mcestimates(nprm, std::vector<double>(N) ), deviances(N) {}
00031 PsiMClist ( const PsiMClist& mclist ) : mcestimates ( mclist.mcestimates ), deviances ( mclist.deviances ) {}
00032 ~PsiMClist ( ) {}
00033 std::vector<double> getEst ( unsigned int i ) const;
00034 double getEst (
00035 unsigned int i,
00036 unsigned int prm
00037 ) const;
00038 void setEst (
00039 unsigned int i,
00040 const std::vector<double> est,
00041 double deviance
00042 );
00043 virtual void setdeviance ( unsigned int i, double deviance );
00044 virtual double getPercentile (
00045 double p,
00046 unsigned int prm
00047 );
00048 virtual double getMean (
00049 unsigned int prm
00050 ) const ;
00051 virtual double getStd (
00052 unsigned int prm
00053 ) const ;
00054 double getdeviance ( unsigned int i ) const;
00055 unsigned int getNsamples ( void ) const { return mcestimates[0].size(); }
00056 unsigned int getNparams ( void ) const { return mcestimates.size(); }
00057 double getDeviancePercentile ( double p );
00058 };
00059
00069 class BootstrapList : public PsiMClist
00070 {
00071 private:
00072 bool BCa;
00073 std::vector<double> acceleration_t;
00074 std::vector<double> bias_t;
00075 std::vector<double> acceleration_s;
00076 std::vector<double> bias_s;
00077 std::vector< std::vector<int> > data;
00078 std::vector<double> cuts;
00079 std::vector< std::vector<double> > thresholds;
00080 std::vector< std::vector<double> > slopes;
00081 std::vector<double> Rpd;
00082 std::vector<double> Rkd;
00083 public:
00084 BootstrapList (
00085 unsigned int N,
00086 unsigned int nprm,
00087 unsigned int nblocks,
00088 std::vector<double> Cuts
00089 ) : PsiMClist (N,nprm),
00090 BCa(false),
00091 acceleration_t(Cuts.size()),
00092 bias_t(Cuts.size()),
00093 acceleration_s(Cuts.size()),
00094 bias_s(Cuts.size()),
00095 data(N,std::vector<int>(nblocks)),
00096 cuts(Cuts),
00097 thresholds (Cuts.size(), std::vector<double> (N)),
00098 slopes (Cuts.size(), std::vector<double> (N)),
00099 Rpd(N),
00100 Rkd(N)
00101 { }
00102
00103 void setBCa_t (
00104 unsigned int i,
00105 double Bias,
00106 double Acceleration
00107 ) { BCa=true; bias_t[i] = Bias; acceleration_t[i] = Acceleration; }
00108 void setBCa_s (
00109 unsigned int i,
00110 double Bias,
00111 double Acceleration
00112 ) { BCa=true; bias_s[i] = Bias; acceleration_s[i] = Acceleration; }
00113 void setData (
00114 unsigned int i,
00115 const std::vector<int>& newdata
00116 );
00117 std::vector<int> getData ( unsigned int i ) const;
00118
00119 double getThres ( double p, unsigned int cut );
00120 double getThres_byPos ( unsigned int i, unsigned int cut );
00121 void setThres ( double thres,
00122 unsigned int i,
00123 unsigned int cut
00124 );
00125
00126 double getSlope ( double p, unsigned int cut );
00127 double getSlope_byPos ( unsigned int i, unsigned int cut );
00128 void setSlope ( double sl,
00129 unsigned int i,
00130 unsigned int cut
00131 );
00132
00133 unsigned int getNblocks ( void ) const { return data[0].size(); }
00134 double getCut ( unsigned int i ) const;
00135 double getAcc_t ( unsigned int i ) const { return acceleration_t[i]; }
00136 double getBias_t ( unsigned int i ) const { return bias_t[i]; }
00137 double getAcc_s ( unsigned int i ) const { return acceleration_s[i]; }
00138 double getBias_s ( unsigned int i ) const { return bias_s[i]; }
00139
00140 void setRpd ( unsigned int i,
00141 double r_pd
00142 );
00143 double getRpd ( unsigned int i ) const;
00144 double percRpd ( double p );
00145
00146 void setRkd ( unsigned int i, double r_kd );
00147 double getRkd ( unsigned int i ) const;
00148 double percRkd ( double p );
00149 };
00150
00156 class JackKnifeList : public PsiMClist
00157 {
00158 private:
00159 double maxdeviance;
00160 std::vector<double> mlestimate;
00161 public:
00162 JackKnifeList (
00163 unsigned int nblocks,
00164 unsigned int nprm,
00165 double maxldev,
00166 std::vector<double> maxlest
00167 ) : PsiMClist ( nblocks, nprm ), maxdeviance(maxldev), mlestimate(maxlest) {}
00168 unsigned int getNblocks ( void ) const { return getNsamples(); }
00169
00180 double influential ( unsigned int block, const std::vector<double>& ci_lower, const std::vector<double>& ci_upper ) const;
00188 bool outlier ( unsigned int block ) const ;
00189 };
00190
00199 class MCMCList : public PsiMClist
00200 {
00201 private:
00202 std::vector<double> posterior_Rpd;
00203 std::vector<double> posterior_Rkd;
00204 std::vector< std::vector<int> > posterior_predictive_data;
00205 std::vector<double> posterior_predictive_deviances;
00206 std::vector<double> posterior_predictive_Rpd;
00207 std::vector<double> posterior_predictive_Rkd;
00208 std::vector< std::vector<double> > logratios;
00209 double accept_rate;
00210 public:
00211 MCMCList (
00212 unsigned int N,
00213 unsigned int nprm,
00214 unsigned int nblocks
00215 ) : PsiMClist ( N, nprm),
00216 posterior_Rpd(N),
00217 posterior_Rkd(N),
00218 posterior_predictive_data(N,std::vector<int>(nblocks)),
00219 posterior_predictive_deviances ( N ),
00220 posterior_predictive_Rpd ( N ),
00221 posterior_predictive_Rkd ( N ),
00222 logratios ( N, std::vector<double>(nblocks ) ) {};
00223 void setppData (
00224 unsigned int i,
00225 const std::vector<int>& ppdata,
00226 double ppdeviance
00227 );
00228 std::vector<int> getppData ( unsigned int i ) const;
00229 int getppData ( unsigned int i, unsigned int j ) const;
00230 double getppDeviance ( unsigned int i ) const;
00231 void setppRpd ( unsigned int i, double Rpd );
00232 double getppRpd ( unsigned int i ) const;
00233 void setppRkd ( unsigned int i, double Rkd );
00234 double getppRkd ( unsigned int i ) const;
00235 void setRpd ( unsigned int i, double Rpd );
00236 double getRpd ( unsigned int i ) const;
00237 void setRkd ( unsigned int i, double Rkd );
00238 double getRkd ( unsigned int i ) const;
00239 unsigned int getNblocks ( void ) const { return posterior_predictive_data[0].size(); }
00240 void setlogratio ( unsigned int i, unsigned int j, double logratio );
00241 double getlogratio ( unsigned int i, unsigned int j ) const;
00242 void set_accept_rate(double rate) {accept_rate = rate; }
00243 double get_accept_rate(void) const {return accept_rate; }
00244 };
00245
00246 void newsample ( const PsiData * data, const std::vector<double>& p, std::vector<int> * sample );
00247
00248 #endif
