SedSat3 1.1.6
Sediment Source Apportionment Tool - Advanced statistical methods for environmental pollution research
Loading...
Searching...
No Matches
cmbdistribution.cpp
Go to the documentation of this file.
1#include "cmbdistribution.h"
2#include "math.h"
3
4double Distribution::pi = 4 * atan(1.0);
5
10
11double Distribution::Eval(const double &x) const
12{
13 double out = 0;
14 if (distribution==distribution_type::normal)
15 out = 1.0 / (sqrt(2*pi)*parameters[1])*exp(-pow(x - parameters[0], 2) / (2 * parameters[1]*parameters[1]));
16 if (distribution==distribution_type::lognormal)
17 out = 1.0 / (sqrt(2*pi)*parameters[1] * x)*exp(-pow(log(x) - parameters[0], 2) / (2 * parameters[1] * parameters[1]));
18 if (distribution==distribution_type::uniform || distribution==distribution_type::dirichlet)
19 if (x<0 || x>1)
20 out = 1e-32;
21 return out;
22}
23
24double Distribution::Eval(const double &x, const vector<double> parameters, distribution_type distribution)
25{
26 double out = 0;
27 if (distribution==distribution_type::normal)
28 out = 1.0 / (sqrt(2*pi)*parameters[1])*exp(-pow(x - parameters[0], 2) / (2 * parameters[1]*parameters[1]));
29 if (distribution==distribution_type::lognormal)
30 out = 1.0 / (sqrt(2*pi)*parameters[1] * x)*exp(-pow(log(x) - parameters[0], 2) / (2 * parameters[1] * parameters[1]));
31 if (distribution==distribution_type::uniform || distribution==distribution_type::dirichlet)
32 if (x<0 || x>1)
33 out = 1e-32;
34 return out;
35}
36
37double Distribution::EvalLog(const double &x)
38{
39 double out = 0;
40 if (distribution==distribution_type::normal)
41 out = -log(sqrt(2*pi)*parameters[1]) -pow(x - parameters[0], 2) / (2 * parameters[1]*parameters[1]);
42 if (distribution==distribution_type::lognormal)
43 out = -log(sqrt(2*pi)*parameters[1]*x) -pow(log(x) - parameters[0], 2) / (2 * parameters[1] * parameters[1]);
44 if (distribution==distribution_type::uniform || distribution==distribution_type::dirichlet)
45 if (x<0 || x>1)
46 out = -1e6;
47
48 return out;
49}
50
51TimeSeries<double> Distribution::EvaluateAsTimeSeries(int numberofpoints, const double &stdcoeff)
52{
53 TimeSeries<double> out;
54 double x;
55 for (int i=0; i<numberofpoints; i++)
56 {
57 if (distribution == distribution_type::normal)
58 {
59 x = parameters[0] - stdcoeff*parameters[1] + i*2*stdcoeff*parameters[1]/(numberofpoints-1);
60 }
61 else if (distribution == distribution_type::lognormal)
62 {
63 x = exp(parameters[0] - stdcoeff*parameters[1] + i*2*stdcoeff*parameters[1]/(numberofpoints-1));
64 }
65 out.append(x,Eval(x));
66 }
67 return out;
68}
69
70TimeSeries<double> Distribution::EvaluateAsTimeSeries(int numberofpoints, const double &stdcoeff, const vector<double> parameters, distribution_type &dist_type)
71{
72 TimeSeries<double> out;
73 double x;
74 for (int i=0; i<numberofpoints; i++)
75 {
76 if (dist_type == distribution_type::normal)
77 {
78 x = parameters[0] - stdcoeff*parameters[1] + i*2*stdcoeff*parameters[1]/(numberofpoints-1);
79 }
80 else if (dist_type == distribution_type::lognormal)
81 {
82 x = exp(parameters[0] - stdcoeff*parameters[1] + i*2*stdcoeff*parameters[1]/(numberofpoints-1));
83 }
84 out.append(x,Eval(x,parameters,dist_type));
85 }
86 return out;
87}
88
89
90Distribution::Distribution(const Distribution &dist)
91{
92 pi = 4 * atan(1.0);
93 parameters = dist.parameters;
94 distribution = dist.distribution;
95 mean_val = dist.mean_val;
96 std_val = dist.std_val;
97 SetType(dist.distribution);
98}
99Distribution& Distribution::operator = (const Distribution &dist)
100{
101 pi = 4 * atan(1.0);
102 parameters = dist.parameters;
103 distribution = dist.distribution;
104 mean_val = dist.mean_val;
105 std_val = dist.std_val;
106 SetType(dist.distribution);
107 return *this;
108}
109
110void Distribution::SetType(const distribution_type &typ)
111{
112 distribution = typ;
113 if (typ == distribution_type::lognormal)
114 parameters.resize(2);
115 if (typ == distribution_type::normal)
116 parameters.resize(2);
117}
118
119
120double Distribution::Mean(parameter_mode param_mode)
121{
122 if (param_mode==parameter_mode::direct)
123 return mean_val;
124
125 if (distribution == distribution_type::lognormal)
126 return exp(parameters[0]+pow(parameters[1],2)/2);
127 if (distribution == distribution_type::normal)
128 return parameters[0];
129
130 return 0;
131}
132
133bool Distribution::operator==(const string &dist_type)
134{
135 if (distribution==distribution_type::normal && dist_type == "normal")
136 return true;
137 if (distribution==distribution_type::lognormal && dist_type == "log-normal")
138 return true;
139 if (distribution==distribution_type::uniform && dist_type == "uniform")
140 return true;
141 if (distribution==distribution_type::dirichlet && dist_type == "uniform")
142 return true;
143 if (distribution==distribution_type::dirichlet && dist_type == "dirichlet")
144 return true;
145
146 return false;
147
148}
149
Distribution
Represents a parametric probability distribution for uncertainty quantification.
Definition cmbdistribution.h:96
Distribution::Eval
double Eval(const double &x) const
Evaluate probability density function at a given value.
Definition cmbdistribution.cpp:11
Distribution::parameters
vector< double > parameters
Distribution parameters vector.
Definition cmbdistribution.h:267
Distribution::EvalLog
double EvalLog(const double &x)
Evaluate natural logarithm of probability density.
Definition cmbdistribution.cpp:37
Distribution::std_val
double std_val
Empirical standard deviation from original sample data.
Definition cmbdistribution.h:399
Distribution::EvaluateAsTimeSeries
TimeSeries< double > EvaluateAsTimeSeries(int numberofpoint=100, const double &stdcoeff=4)
Generate time series representation of the distribution for plotting.
Definition cmbdistribution.cpp:51
Distribution::distribution
distribution_type distribution
The type of probability distribution.
Definition cmbdistribution.h:275
Distribution::operator=
Distribution & operator=(const Distribution &dist)
Assignment operator.
Definition cmbdistribution.cpp:99
Distribution::SetType
void SetType(const distribution_type &typ)
Set the distribution type and resize parameters vector.
Definition cmbdistribution.cpp:110
Distribution::pi
static double pi
Mathematical constant π (pi)
Definition cmbdistribution.h:380
Distribution::mean_val
double mean_val
Empirical mean from original sample data.
Definition cmbdistribution.h:391
Distribution::Mean
double Mean(parameter_mode param_mode=parameter_mode::based_on_fitted_distribution)
Calculate the mean (expected value) of the distribution.
Definition cmbdistribution.cpp:120
Distribution::operator==
bool operator==(const string &dist_type)
Compare distribution type with a string identifier.
Definition cmbdistribution.cpp:133
Distribution::Distribution
Distribution()
Default constructor.
Definition cmbdistribution.cpp:6
parameter_mode
parameter_mode
Specifies whether to use direct data statistics or fitted distribution parameters.
Definition cmbdistribution.h:30
parameter_mode::direct
@ direct
Use empirical statistics from data.
distribution_type
distribution_type
Enumeration of probability distribution types supported in SedSat3.
Definition cmbdistribution.h:14
distribution_type::dirichlet
@ dirichlet
Dirichlet distribution (treated as uniform on simplex in current implementation)
distribution_type::lognormal
@ lognormal
Lognormal distribution: ln(x) ~ N(μ, σ²), for strictly positive variables.
distribution_type::uniform
@ uniform
Uniform distribution on [0,1].
distribution_type::normal
@ normal
Normal (Gaussian) distribution: p(x) = N(μ, σ²)