1+ #ifndef PROBABILITY_H
2+ #define PROBABILITY_H
3+ #include " numerical.h"
4+ using namespace nmr ;
5+
6+ namespace prb {
7+
8+ // -----声明部分-----
9+
10+ // 随机变量
11+ class DiscreteRandomVariable
12+ {
13+ public:
14+ DiscreteRandomVariable (const Matrix<double >&);
15+ friend Matrix<double > Distribution (const DiscreteRandomVariable &);
16+ private:
17+ Matrix<double > _distribution;
18+ };
19+ Matrix<double > Distribution (const DiscreteRandomVariable &);
20+ double Expectation (Matrix<double > );
21+ double Expectation (DiscreteRandomVariable X);
22+
23+ // 分布
24+ class NormalDistribution
25+ {
26+ public:
27+ NormalDistribution (double , double );
28+ friend double Expectation (const NormalDistribution&);
29+ friend double StandardDeviation (const NormalDistribution&);
30+
31+ private:
32+ double _mu;
33+ double _sigma;
34+ };
35+ double Expectation (const NormalDistribution &);
36+ double StandardDeviation (const NormalDistribution &);
37+ double PDF (const NormalDistribution &, double );
38+ double CDF (const NormalDistribution &, double );
39+
40+ // -----定义部分-----
41+
42+ // 随机变量
43+ DiscreteRandomVariable::DiscreteRandomVariable (const Matrix<double > &M)
44+ {
45+ _distribution = M;
46+ }
47+ Matrix<double > Distribution (const DiscreteRandomVariable &X)
48+ {
49+ return X._distribution ;
50+ }
51+ double Expectation (Matrix<double > A)
52+ {
53+ if (RowSize (A)!=2 )
54+ {
55+ cerr << " 错误:输入的矩阵不是分布列" << ' \n ' ;
56+ return 0 ;
57+ }
58+ double s = 0 ;
59+ for (int j= 0 ; j< ColumnSize (A);++j)
60+ {
61+ s = s + A[0 ][j] * A[1 ][j];
62+ }
63+ return s;
64+ }
65+ double Expectation (DiscreteRandomVariable X)
66+ {
67+ return Expectation (Distribution (X));
68+ }
69+
70+ // 分布
71+ NormalDistribution::NormalDistribution (double mu,double sigma)
72+ {
73+ _mu = mu;
74+ _sigma = sigma;
75+ }
76+ double Expectation (const NormalDistribution &N)
77+ {
78+ return N._mu ;
79+ }
80+ double StandardDeviation (const NormalDistribution &N)
81+ {
82+ return N._sigma ;
83+ }
84+ double PDF (const NormalDistribution &N,double x)
85+ {
86+ const double mu = Expectation (N);
87+ const double sigma = StandardDeviation (N);
88+ return 1 / (sqrt (2 * Pi<double >)*sigma)*exp (-pow (x-mu,2 )/(2 *sigma*sigma));
89+ }
90+ double CDF (const NormalDistribution &N,double x)
91+ {
92+ double leftend = Expectation (N) - 6 * StandardDeviation (N);
93+ if (x<=leftend)
94+ return 0 ;
95+ auto p = [&](double t) { return PDF (N, t); };
96+ return Integrate<double >(p, leftend, x);
97+ }
98+
99+ }
100+
101+ #endif
0 commit comments