1 /++ 2 This module contains algorithms for the $(LINK2 https://en.wikipedia.org/wiki/Pareto_distribution, Pareto Distribution). 3 4 License: $(HTTP www.apache.org/licenses/LICENSE-2.0, Apache-2.0) 5 6 Authors: John Michael Hall 7 8 Copyright: 2022-3 Mir Stat Authors. 9 10 +/ 11 12 module mir.stat.distribution.pareto; 13 14 import mir.internal.utility: isFloatingPoint; 15 16 /++ 17 Computes the pareto probability density function (PDF). 18 19 Params: 20 x = value to evaluate PDF 21 xMin = scale parameter 22 alpha = shape parameter 23 24 See_also: 25 $(LINK2 https://en.wikipedia.org/wiki/Pareto_distribution, Pareto Distribution) 26 +/ 27 @safe pure nothrow @nogc 28 T paretoPDF(T)(const T x, const T xMin, const T alpha) 29 if (isFloatingPoint!T) 30 in (x >= xMin, "x must be greater than or equal to xMin") 31 in (xMin > 0, "xMin must be greater than zero") 32 in (alpha > 0, "alpha must be greater than zero") 33 { 34 import mir.math.common: pow; 35 36 return alpha * pow(xMin, alpha) / pow(x, alpha + 1); 37 } 38 39 /// 40 version(mir_stat_test) 41 @safe pure nothrow @nogc 42 unittest { 43 import mir.test: shouldApprox; 44 45 1.0.paretoPDF(1, 3).shouldApprox == 3; 46 2.0.paretoPDF(1, 3).shouldApprox == 0.1875; 47 3.0.paretoPDF(2, 4).shouldApprox == 0.2633745; 48 } 49 50 /++ 51 Computes the pareto cumulative distribution function (CDF). 52 53 Params: 54 x = value to evaluate CDF 55 xMin = scale parameter 56 alpha = shape parameter 57 58 See_also: 59 $(LINK2 https://en.wikipedia.org/wiki/Pareto_distribution, Pareto Distribution) 60 +/ 61 @safe pure nothrow @nogc 62 T paretoCDF(T)(const T x, const T xMin, const T alpha) 63 if (isFloatingPoint!T) 64 in (x >= xMin, "x must be greater than or equal to xMin") 65 in (xMin > 0, "xMin must be greater than zero") 66 in (alpha > 0, "alpha must be greater than zero") 67 { 68 import mir.math.common: pow; 69 70 return 1 - pow(xMin / x, alpha); 71 } 72 73 /// 74 version(mir_stat_test) 75 @safe pure nothrow @nogc 76 unittest { 77 import mir.test: shouldApprox; 78 79 1.0.paretoCDF(1, 3).shouldApprox == 0; 80 2.0.paretoCDF(1, 3).shouldApprox == 0.875; 81 3.0.paretoCDF(2, 4).shouldApprox == 0.8024691; 82 } 83 84 /++ 85 Computes the pareto complementary cumulative distribution function (CCDF). 86 87 Params: 88 x = value to evaluate CCDF 89 xMin = scale parameter 90 alpha = shape parameter 91 92 See_also: 93 $(LINK2 https://en.wikipedia.org/wiki/Pareto_distribution, Pareto Distribution) 94 +/ 95 @safe pure nothrow @nogc 96 T paretoCCDF(T)(const T x, const T xMin, const T alpha) 97 if (isFloatingPoint!T) 98 in (x >= xMin, "x must be greater than or equal to xMin") 99 in (xMin > 0, "xMin must be greater than zero") 100 in (alpha > 0, "alpha must be greater than zero") 101 { 102 import mir.math.common: pow; 103 104 return pow(xMin / x, alpha); 105 } 106 107 /// 108 version(mir_stat_test) 109 @safe pure nothrow @nogc 110 unittest { 111 import mir.test: shouldApprox; 112 113 1.0.paretoCCDF(1, 3).shouldApprox == 1; 114 2.0.paretoCCDF(1, 3).shouldApprox == 0.125; 115 3.0.paretoCCDF(2, 4).shouldApprox == 0.1975309; 116 } 117 118 /++ 119 Computes the pareto inverse cumulative distribution function (InvCDF). 120 121 Params: 122 p = value to evaluate InvCDF 123 xMin = scale parameter 124 alpha = shape parameter 125 126 See_also: 127 $(LINK2 https://en.wikipedia.org/wiki/Pareto_distribution, Pareto Distribution) 128 +/ 129 @safe pure nothrow @nogc 130 T paretoInvCDF(T)(const T p, const T xMin, const T alpha) 131 if (isFloatingPoint!T) 132 in (p >= 0, "p must be greater than or equal to 0") 133 in (p <= 1, "p must be less than or equal to 1") 134 in (xMin > 0, "xMin must be greater than zero") 135 in (alpha > 0, "alpha must be greater than zero") 136 { 137 import mir.math.common: pow; 138 139 return xMin / pow(1 - p, cast(T) 1 / alpha); 140 } 141 142 /// 143 version(mir_stat_test) 144 @safe pure nothrow @nogc 145 unittest { 146 import mir.test: shouldApprox; 147 148 0.0.paretoInvCDF(1, 3).shouldApprox == 1; 149 0.875.paretoInvCDF(1, 3).shouldApprox == 2; 150 0.8024691.paretoInvCDF(2, 4).shouldApprox == 3; 151 } 152 153 /++ 154 Computes the pareto log probability density function (LPDF). 155 156 Params: 157 x = value to evaluate LPDF 158 xMin = scale parameter 159 alpha = shape parameter 160 161 See_also: 162 $(LINK2 https://en.wikipedia.org/wiki/Pareto_distribution, Pareto Distribution) 163 +/ 164 @safe pure nothrow @nogc 165 T paretoLPDF(T)(const T x, const T xMin, const T alpha) 166 if (isFloatingPoint!T) 167 in (x >= xMin, "x must be greater than or equal to xMin") 168 in (xMin > 0, "xMin must be greater than zero") 169 in (alpha > 0, "alpha must be greater than zero") 170 { 171 import mir.math.common: log; 172 import mir.math.internal.xlogy: xlogy; 173 174 return log(alpha) + xlogy(alpha, xMin) - xlogy(alpha + 1, x); 175 } 176 177 /// 178 version(mir_stat_test) 179 @safe pure nothrow @nogc 180 unittest { 181 import mir.math.common: log; 182 import mir.test: shouldApprox; 183 184 1.0.paretoLPDF(1, 3).shouldApprox == log(paretoPDF(1.0, 1, 3)); 185 2.0.paretoLPDF(1, 3).shouldApprox == log(paretoPDF(2.0, 1, 3)); 186 3.0.paretoLPDF(2, 4).shouldApprox == log(paretoPDF(3.0, 2, 4)); 187 }