root/library/tests/epdf_harness.cpp @ 507

Revision 507, 7.9 kB (checked in by vbarta, 15 years ago)

removed class compositepdf; keeping mpdfs of mprod and merger_base in shared pointers

Line 
1#include "epdf_harness.h"
2#include "base/bdmbase.h"
3#include "base/user_info.h"
4#include "stat/exp_family.h"
5#include "mat_checks.h"
6#include "test_util.h"
7#include "UnitTest++.h"
8#include <memory>
9
10namespace bdm {
11
12void epdf_harness::test_config ( const char *config_file_name ) {
13        RV::clear_all();
14
15        UIFile in ( config_file_name );
16        Array<epdf_harness *> input;
17        UI::get ( input, in, "data", UI::compulsory );
18        int sz = input.size();
19        CHECK ( sz > 0 );
20        for ( int i = 0; i < sz; ++i ) {
21                input ( i )->test ( config_file_name, i );
22        }
23}
24
25void epdf_harness::from_setting ( const Setting &set ) {
26        hepdf = UI::build<epdf> ( set, "epdf", UI::compulsory );
27        UI::get ( mean, set, "mean", UI::compulsory );
28        UI::get ( variance, set, "variance", UI::compulsory );
29
30        UI::get (support, set, "support", UI::optional );
31        UI::get (nbins, set, "nbins", UI::optional );
32        UI::get (nsamples, set, "nsamples", UI::optional );
33        UI::get (R, set, "R", UI::optional );
34
35        mrv = UI::build<RV> (set, "marginal_rv", UI::optional );
36
37        UI::get ( tolerance, set, "tolerance", UI::optional );
38}
39
40void epdf_harness::test ( const char *config_name, int idx ) {
41        CurrentContext cc ( config_name, idx );
42
43        CHECK_CLOSE_EX ( mean, hepdf->mean(), tolerance );
44        CHECK_CLOSE_EX ( variance, hepdf->variance(), tolerance );
45
46        if ( support.rows() == 2 ) {
47                int old_size = nbins.size();
48                if ( old_size < 2 ) {
49                        ivec new_nbins ( "100 100" );
50                        for ( int i = 0; i < old_size; ++i ) {
51                                new_nbins ( i ) = nbins ( i );
52                        }
53
54                        nbins = new_nbins;
55                }
56
57                check_support_mean();
58                check_support_integral();
59        }
60
61        if ( R.rows() > 0 ) {
62                check_sample_mean();
63                check_covariance();
64        }
65
66        if ( mrv.get() ) {
67                RV crv = hepdf->_rv().subt ( *mrv );
68                shared_ptr<epdf> m = hepdf->marginal ( *mrv );
69                shared_ptr<mpdf> c = hepdf->condition ( crv );
70
71                Array<shared_ptr<mpdf> > aa ( 2 );
72                aa ( 0 ) = c;
73                aa ( 1 ) = new mepdf ( m );
74                mprod mEp ( aa );
75
76                check_cond_mean(mEp);
77                if ( R.rows() > 0 ) {
78                        check_cond_covariance( mEp );
79                }
80
81                // test of pdflog at zero
82                vec zero ( 0 );
83                vec zeron ( hepdf->dimension() );
84                for ( int i = 0; i < zeron.size(); ++i ) {
85                        zeron ( i ) = 0;
86                }
87
88                double lpz = hepdf->evallog ( zeron );
89                CHECK_CLOSE_EX ( lpz, mEp.evallogcond ( zeron, zero ), tolerance );
90
91                vec lpzv(1);
92                lpzv(0) = lpz;
93
94                mat zero1n ( hepdf->dimension(), 1 );
95                for ( int i = 0; i < zero1n.rows(); ++i ) {
96                        zero1n ( i, 0 ) = 0;
97                }
98
99                vec lpzv_act = hepdf->evallog_m ( zero1n );
100                CHECK_CLOSE_EX ( lpzv, lpzv_act, tolerance );
101
102                Array<vec> zeroa(3);
103                lpzv = vec( zeroa.size() );
104                for ( int i = 0; i < zeroa.size(); ++i ) {
105                        zeroa(i) = zeron;
106                        lpzv(i) = lpz;
107                }
108
109                lpzv_act = hepdf->evallog_m ( zeroa );
110                CHECK_CLOSE_EX ( lpzv, lpzv_act, tolerance );
111        }
112}
113
114void epdf_harness::check_support_mean() {
115        vec xb = support.get_row ( 0 );
116        vec yb = support.get_row ( 1 );
117
118        int tc = 0;
119        Array<vec> actual(CurrentContext::max_trial_count);
120        do {
121                vec emu = num_mean2 ( hepdf.get(), xb, yb, nbins ( 0 ), nbins ( 1 ) );
122                actual( tc ) = emu;
123                ++tc;
124        } while ( ( tc < CurrentContext::max_trial_count ) &&
125                  !UnitTest::AreClose ( mean, actual( tc - 1 ), tolerance ) );
126        if ( ( tc == CurrentContext::max_trial_count ) &&
127             ( !UnitTest::AreClose ( mean, actual( CurrentContext::max_trial_count - 1 ), tolerance ) ) ) {
128                UnitTest::MemoryOutStream stream;
129                stream << CurrentContext::format_context(__LINE__) << "expected " << mean << " +/- " << tolerance << " but was " << actual;
130
131                UnitTest::TestDetails details(*UnitTest::CurrentTest::Details(), 0, false);
132
133                UnitTest::CurrentTest::Results()->OnTestFailure ( details, stream.GetText() );
134       }
135}
136
137void epdf_harness::check_support_integral() {
138        vec xb = support.get_row ( 0 );
139        vec yb = support.get_row ( 1 );
140
141        int tc = 0;
142        Array<double> actual(CurrentContext::max_trial_count);
143        do {
144                double nc = normcoef ( hepdf.get(), xb, yb, nbins ( 0 ), nbins ( 1 ) );
145                actual( tc ) = nc;
146                ++tc;
147        } while ( ( tc < CurrentContext::max_trial_count ) &&
148                  !UnitTest::AreClose ( 1.0, actual( tc - 1 ), tolerance ) );
149        if ( ( tc == CurrentContext::max_trial_count ) &&
150             ( !UnitTest::AreClose ( 1.0, actual( CurrentContext::max_trial_count - 1 ), tolerance ) ) ) {
151                UnitTest::MemoryOutStream stream;
152                stream << CurrentContext::format_context(__LINE__) << "expected " << mean << " +/- " << tolerance << " but was " << actual;
153
154                UnitTest::TestDetails details(*UnitTest::CurrentTest::Details(), 0, false);
155
156                UnitTest::CurrentTest::Results()->OnTestFailure ( details, stream.GetText() );
157       }
158}
159
160void epdf_harness::check_sample_mean() {
161        vec delta = make_close_tolerance ( variance, nsamples );
162
163        int tc = 0;
164        Array<vec> actual(CurrentContext::max_trial_count);
165        do {
166                mat smp = hepdf->sample_m ( nsamples );
167                vec emu = smp * ones ( nsamples ) / nsamples;
168                actual( tc ) = emu;
169                ++tc;
170        } while ( ( tc < CurrentContext::max_trial_count ) &&
171                  !UnitTest::AreClose ( mean, actual( tc - 1 ), delta ) );
172        if ( ( tc == CurrentContext::max_trial_count ) &&
173             ( !UnitTest::AreClose ( mean, actual( CurrentContext::max_trial_count - 1 ), delta ) ) ) {
174                UnitTest::MemoryOutStream stream;
175                stream << CurrentContext::format_context(__LINE__) << "expected " << mean << " +/- " << delta << " but was " << actual;
176
177                UnitTest::TestDetails details(*UnitTest::CurrentTest::Details(), 0, false);
178
179                UnitTest::CurrentTest::Results()->OnTestFailure ( details, stream.GetText() );
180       }
181}
182
183void epdf_harness::check_covariance() {
184        int tc = 0;
185        Array<mat> actual(CurrentContext::max_trial_count);
186        do {
187                mat smp = hepdf->sample_m ( nsamples );
188                vec emu = smp * ones ( nsamples ) / nsamples;
189                mat er = ( smp * smp.T() ) / nsamples - outer_product ( emu, emu );
190                actual( tc ) = er;
191                ++tc;
192        } while ( ( tc < CurrentContext::max_trial_count ) &&
193                  !UnitTest::AreClose ( R, actual( tc - 1 ), tolerance ) );
194       
195        if ( ( tc == CurrentContext::max_trial_count ) &&
196             ( !UnitTest::AreClose ( R, actual( CurrentContext::max_trial_count - 1 ), tolerance ) ) ) {
197                UnitTest::MemoryOutStream stream;
198                stream << CurrentContext::format_context(__LINE__) << "expected " << R << " +/- " << tolerance << " but was " << actual;
199
200                UnitTest::TestDetails details(*UnitTest::CurrentTest::Details(), 0, false);
201
202                UnitTest::CurrentTest::Results()->OnTestFailure ( details, stream.GetText() );
203       }
204}
205
206void epdf_harness::check_cond_mean( mprod &mep ) {
207        vec delta = make_close_tolerance ( variance, nsamples );
208
209        int tc = 0;
210        Array<vec> actual(CurrentContext::max_trial_count);
211        do {
212                mat smp = mep.samplecond ( vec ( 0 ), nsamples );
213                vec emu = sum ( smp, 2 ) / nsamples;
214                actual( tc ) = emu;
215                ++tc;
216        } while ( ( tc < CurrentContext::max_trial_count ) &&
217                  !UnitTest::AreClose ( mean, actual( tc - 1 ), delta ) );
218        if ( ( tc == CurrentContext::max_trial_count ) &&
219             ( !UnitTest::AreClose ( mean, actual( CurrentContext::max_trial_count - 1 ), delta ) ) ) {
220                UnitTest::MemoryOutStream stream;
221                stream << CurrentContext::format_context(__LINE__) << "expected " << mean << " +/- " << delta << " but was " << actual;
222
223                UnitTest::TestDetails details(*UnitTest::CurrentTest::Details(), 0, false);
224
225                UnitTest::CurrentTest::Results()->OnTestFailure ( details, stream.GetText() );
226       }
227}
228
229void epdf_harness::check_cond_covariance( mprod &mep ) {
230        int tc = 0;
231        Array<mat> actual(CurrentContext::max_trial_count);
232        do {
233                mat smp = mep.samplecond ( vec ( 0 ), nsamples );
234                vec emu = sum ( smp, 2 ) / nsamples;
235                mat er = ( smp * smp.T() ) / nsamples - outer_product ( emu, emu );
236                actual( tc ) = er;
237                ++tc;
238        } while ( ( tc < CurrentContext::max_trial_count ) &&
239                  !UnitTest::AreClose ( R, actual( tc - 1 ), tolerance ) );
240        if ( ( tc == CurrentContext::max_trial_count ) &&
241             ( !UnitTest::AreClose ( R, actual( CurrentContext::max_trial_count - 1 ), tolerance ) ) ) {
242                UnitTest::MemoryOutStream stream;
243                stream << CurrentContext::format_context(__LINE__) << "expected " << mean << " +/- " << tolerance << " but was " << actual;
244
245                UnitTest::TestDetails details(*UnitTest::CurrentTest::Details(), 0, false);
246
247                UnitTest::CurrentTest::Results()->OnTestFailure ( details, stream.GetText() );
248       }
249}
250
251}
Note: See TracBrowser for help on using the browser.