00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018 #include "Scaling.h"
00019 #include "TargetInfo.h"
00020
00021 using namespace std;
00022
00023 Scaling slope(const std::valarray<double>& x, const TargetInfo& targets) {
00024
00025 double xx = 0.0;
00026 double xy = 0.0;
00027
00028 const valarray<double>& y = targets.targets();
00029
00030 for (unsigned i = 0; i < x.size(); ++i) {
00031 xx += x[i] * x[i];
00032 xy += x[i] * y[i];
00033 }
00034
00035 if (xx < 1e-7) return Scaling(new LinearScaling(0.0,0.0));
00036
00037 double b = xy / xx;
00038
00039 return Scaling(new LinearScaling(0.0, b));
00040
00041 }
00042
00043
00044 Scaling regularized_least_squares(const std::valarray<double>& inputs, const TargetInfo& targets, double lambda) {
00045
00046 double n = inputs.size();
00047
00048 valarray<double> x = inputs;
00049
00050 double a,b,d;
00051 a=b=d=0;
00052
00053 for (unsigned i = 0; i < n; ++i) {
00054 a += 1 + lambda;
00055 b += x[i];
00056 d += x[i] * x[i] + lambda;
00057 }
00058
00059
00060
00061 double ad_bc = a*d - b * b;
00062
00063
00064 if (ad_bc < 1e-17) return Scaling(new LinearScaling);
00065
00066 double ai = d/ad_bc;
00067 double bi = -b/ad_bc;
00068 double di = a/ad_bc;
00069 double ci = bi;
00070
00071
00072
00073 std::valarray<double> ones = x;
00074
00075
00076 for (unsigned i = 0; i < n; ++i)
00077 {
00078 ones[i] = (ai + bi * x[i]);
00079 x[i] = (ci + di * x[i]);
00080 }
00081
00082
00083
00084 a = 0.0;
00085 b = 0.0;
00086
00087 const valarray<double>& t = targets.targets();
00088
00089 for (unsigned i = 0; i < n; ++i)
00090 {
00091 a += ones[i] * t[i];
00092 b += x[i] * t[i];
00093 }
00094
00095 return Scaling(new LinearScaling(a,b));
00096 }
00097
00098 Scaling ols(const std::valarray<double>& y, const std::valarray<double>& t) {
00099 double n = y.size();
00100
00101 double y_mean = y.sum() / n;
00102 double t_mean = t.sum() / n;
00103
00104 std::valarray<double> y_var = (y - y_mean);
00105 std::valarray<double> t_var = (t - t_mean);
00106 std::valarray<double> cov = t_var * y_var;
00107
00108 y_var *= y_var;
00109 t_var *= t_var;
00110
00111 double sumvar = y_var.sum();
00112
00113 if (sumvar == 0. || sumvar/n < 1e-7 || sumvar/n > 1e+7)
00114 return Scaling(new LinearScaling(t_mean,0.));
00115
00116
00117 double b = cov.sum() / sumvar;
00118 double a = t_mean - b * y_mean;
00119
00120 Scaling s = Scaling(new LinearScaling(a,b));
00121
00122 return s;
00123 }
00124
00125 Scaling ols(const std::valarray<double>& y, const TargetInfo& targets) {
00126 double n = y.size();
00127
00128 double y_mean = y.sum() / n;
00129
00130 std::valarray<double> y_var = (y - y_mean);
00131 std::valarray<double> cov = targets.tcov_part() * y_var;
00132
00133 y_var *= y_var;
00134
00135 double sumvar = y_var.sum();
00136
00137 if (sumvar == 0. || sumvar/n < 1e-7 || sumvar/n > 1e+7)
00138 return Scaling(new LinearScaling(targets.tmean(),0.));
00139
00140
00141 double b = cov.sum() / sumvar;
00142 double a = targets.tmean() - b * y_mean;
00143
00144 if (!finite(b)) {
00145
00146 cout << a << ' ' << b << endl;
00147 cout << sumvar << endl;
00148 cout << y_mean << endl;
00149 cout << cov.sum() << endl;
00150 exit(1);
00151 }
00152
00153 Scaling s = Scaling(new LinearScaling(a,b));
00154
00155 return s;
00156 }
00157
00158
00159 Scaling wls(const std::valarray<double>& inputs, const TargetInfo& targets) {
00160
00161 std::valarray<double> x = inputs;
00162 const std::valarray<double>& w = targets.weights();
00163
00164 unsigned n = x.size();
00165
00166 std::valarray<double> wx = targets.weights() * x;
00167
00168
00169 double a,b,d;
00170 a=b=d=0.0;
00171
00172 for (unsigned i = 0; i < n; ++i)
00173 {
00174 a += w[i];
00175 b += wx[i];
00176 d += x[i] * wx[i];
00177 }
00178
00179
00180
00181 double ad_bc = a*d - b * b;
00182
00183
00184 if (ad_bc < 1e-17) return Scaling(new LinearScaling);
00185
00186 double ai = d/ad_bc;
00187 double bi = -b/ad_bc;
00188 double di = a/ad_bc;
00189 double ci = bi;
00190
00191
00192
00193
00194 std::valarray<double>& ones = wx;
00195
00196
00197 for (unsigned i = 0; i < n; ++i)
00198 {
00199 ones[i] = w[i]*(ai + bi * x[i]);
00200 x[i] = w[i]*(ci + di * x[i]);
00201 }
00202
00203
00204
00205 a = 0.0;
00206 b = 0.0;
00207
00208 const valarray<double>& t = targets.targets();
00209
00210 for (unsigned i = 0; i < n; ++i)
00211 {
00212 a += ones[i] * t[i];
00213 b += x[i] * t[i];
00214 }
00215
00216 return Scaling(new LinearScaling(a,b));
00217 }
00218
00219
00220
00221
00222 double mse(const std::valarray<double>& y, const TargetInfo& t) {
00223
00224 valarray<double> residuals = t.targets()-y;
00225 residuals *= residuals;
00226 double sz = residuals.size();
00227 if (t.has_weights()) {
00228 residuals *= t.weights();
00229 sz = 1.0;
00230 }
00231
00232 return residuals.sum() / sz;
00233 }
00234
00235 double rms(const std::valarray<double>& y, const TargetInfo& t) {
00236 return sqrt(mse(y,t));
00237 }
00238
00239 double mae(const std::valarray<double>& y, const TargetInfo& t) {
00240 valarray<double> residuals = abs(t.targets()-y);
00241 if (t.has_weights()) residuals *= t.weights();
00242 return residuals.sum() / residuals.size();
00243 }
00244
00245
00246
00247
00248
00249
00250
00251
00252
00253
00254
00255
00256
00257
00258
00259
00260
00261
00262
00263
00264
00265
00266
00267
00268
00269
00270
00271
00272
00273
00274
00275
00276
00277
00278
00279
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294
00295
00296
00297
00298
00299
00300
00301
00302
00303
00304
00305
00306
00307
00308
00309
00310
00311
00312
00313
00314
00315
00316
00317
00318
00319
00320
00321
00322
00323
00324
00325
00326
00327
00328
00329
00330
00331
00332
00333
00334
00335
00336
00337
00338
00339
00340
00341
00342
00343
00344
00345
00346
00347
00348
00349
00350
00351
00352
00353
00354
00355
00356
00357
00358
00359
00360
00361
00362
00363
00364
00365
00366
00367
00368
00369
00370
00371
00372
00373
00374
00375
00376
00377
00378
00379
00380
00381
00382
00383
00384
00385
00386
00387
00388
00389
00390
00391
00392
00393
00394
00395
00396
00397
00398
00399
00400
00401
00402
00403
00404
00405
00406
00407
00408
00409
00410
00411
00412
00413
00414
00415
00416
00417