Line data Source code
1 : /**
2 : Bojan Nikolic <bojan@bnikolic.co.uk>
3 : Initial version 2009
4 :
5 : This file is part of BNMin1 and is licensed under GNU General
6 : Public License version 2
7 :
8 : \file mcpoint.hxx
9 : Renamed to mcpoint.h
10 :
11 : Structure representing a Markov Chain point
12 : */
13 :
14 : #ifndef _BNMIN1_MCPOINT_HXX__
15 : #define _BNMIN1_MCPOINT_HXX__
16 :
17 : #include <vector>
18 : #include <list>
19 : #include <set>
20 :
21 : namespace Minim {
22 :
23 : /** \brief Data to be recorded at each point in an MCMC distribution
24 : */
25 : struct MCPoint
26 : {
27 : /// The actual parameters
28 : std::vector<double> p;
29 : /// Log-likelihood of this point
30 : double ll;
31 : /// A vector to store derived quantities at sample of the
32 : /// distribution
33 : std::vector<double> fval;
34 :
35 : /// Default constructor allowed, fill in the data later
36 : MCPoint(void);
37 :
38 : /** \Construct with supplied position vector
39 : */
40 : MCPoint(const std::vector<double> &p);
41 :
42 : /** \brief The parameter vector has n values
43 : */
44 : MCPoint(size_t np);
45 :
46 : MCPoint(const MCPoint &other);
47 :
48 : MCPoint & operator=(const MCPoint &other);
49 :
50 :
51 : };
52 :
53 : /** \brief Define ordering of MCPoints on basis of their likelihood
54 :
55 : Should be useful for nested sampling
56 : */
57 7853994 : inline bool operator< (const MCPoint &a, const MCPoint &b)
58 : {
59 7853994 : return a.ll < b.ll;
60 : }
61 :
62 : /** \brief Weighted posterior point
63 :
64 : As produced by non markov chain posterior exploration
65 : */
66 : struct WPPoint:
67 : public MCPoint
68 : {
69 : /** \brief Weight factor
70 : */
71 : double w;
72 :
73 : WPPoint(void):
74 : w(0.0)
75 : {
76 : }
77 :
78 : WPPoint(const std::vector<double> &p,
79 : double w):
80 : MCPoint(p),
81 : w(w)
82 : {
83 : }
84 :
85 : /** \brief Construct from MCPoint and a supplied weight
86 : */
87 771020 : WPPoint(const MCPoint &mp,
88 771020 : double w):
89 : MCPoint(mp),
90 771020 : w(w)
91 : {
92 771020 : }
93 :
94 : };
95 :
96 : /** \brief Calculate the first moment of each parameter
97 : */
98 : void moment1(const std::list<WPPoint> &l,
99 : std::vector<double> &res);
100 :
101 : void moment1(const std::list<WPPoint> &l,
102 : double Z,
103 : std::vector<double> &res);
104 :
105 : /** \brief Calculate the second moment of each parameter
106 : */
107 : void moment2(const std::list<WPPoint> &l,
108 : const std::vector<double> &m1,
109 : std::vector<double> &res);
110 :
111 : void moment2(const std::list<WPPoint> &l,
112 : const std::vector<double> &m1,
113 : double Z,
114 : std::vector<double> &res);
115 :
116 :
117 : void moment1(const std::set<MCPoint> &s,
118 : std::vector<double> &res);
119 :
120 : void moment2(const std::set<MCPoint> &s,
121 : const std::vector<double> &m1,
122 : std::vector<double> &res);
123 :
124 : /** \brief Outer moment 2, i.e., covariances
125 : */
126 : void omoment2(const std::set<MCPoint> &s,
127 : const std::vector<double> &m1,
128 : std::vector<double> &res);
129 :
130 : /** Covariances, leaving the m1 to be calculated internally
131 : */
132 : void omoment2(const std::set<MCPoint> &s,
133 : std::vector<double> &res);
134 :
135 : /** Standard eviation of each coordiante separately
136 : */
137 : void StdDev(const std::set<MCPoint> &s,
138 : std::vector<double> &res);
139 :
140 :
141 : /** \brief Compute the principal directions from the covariance
142 : matrix
143 :
144 : \param eigvals The eigen values
145 : \param eigvects The eigen vectors
146 :
147 : */
148 : void principalCV(const std::vector<double> &cv,
149 : std::vector<double> &eigvals,
150 : std::vector<double> &eigvects);
151 :
152 :
153 : /** \brief Compute the histogram of the posterior from weighted
154 : sample points
155 :
156 : \param low The low boundary of the region to be binned in each
157 : dimension
158 :
159 : \param high The high boundary of the region to be binned
160 :
161 : \param nbins Number of bins in each direction
162 :
163 : \param res Results are stored here
164 :
165 : */
166 : void postHist(const std::list<WPPoint> &l,
167 : double Z,
168 : const std::vector<double> &low,
169 : const std::vector<double> &high,
170 : size_t nbins,
171 : std::vector<double> &res);
172 :
173 : /** Calculate the marginalised probability histogram
174 :
175 : \param pi Parameter index to marginaise to
176 :
177 : \param Z Evidence value, used to normalize the histogram. If
178 : this is set to one, the integral will be the evidence rather
179 : than the marginalised probability given the hypothesis
180 :
181 : \param high The high boundary of region to histogram, i.e., the
182 : *top* of the highest bin
183 :
184 : */
185 : void marginHist(const std::list<WPPoint> &l,
186 : size_t pi,
187 : double Z,
188 : double low,
189 : double high,
190 : size_t nbins,
191 : std::vector<double> &res);
192 :
193 : /** Calculate the marginalisation of the posterior to two dimensions
194 :
195 : \param Z Evidence value, used to normalize the histogram. If
196 : this is set to one, the integral will be the evidence rather
197 : than the marginalised probability given the hypothesis
198 :
199 : */
200 : void marginHist2D(const std::list<WPPoint> &l,
201 : double Z,
202 : size_t i,
203 : double ilow,
204 : double ihigh,
205 : size_t j,
206 : double jlow,
207 : double jhigh,
208 : size_t nbins,
209 : std::vector<double> &res);
210 :
211 :
212 : }
213 :
214 : #endif
|
