1 module rmathd.lnorm; 2 3 public import rmathd.common; 4 public import rmathd.normal; 5 6 /* 7 ** dmd lnorm.d common.d normal.d && ./lnorm 8 */ 9 10 T dlnorm(T: double)(T x, T meanlog, T sdlog, int give_log) 11 { 12 T y; 13 mixin R_D__0!give_log; 14 15 if (isNaN(x) || isNaN(meanlog) || isNaN(sdlog)) 16 return x + meanlog + sdlog; 17 18 if(sdlog <= 0) { 19 if(sdlog < 0) 20 return T.nan; 21 // sdlog == 0 : 22 return (log(x) == meanlog) ? T.infinity : R_D__0; 23 } 24 if(x <= 0) 25 return R_D__0; 26 27 y = (log(x) - meanlog) / sdlog; 28 return (give_log ? -(M_LN_SQRT_2PI + 0.5 * y * y + log(x * sdlog)) : 29 M_1_SQRT_2PI * exp(-0.5 * y * y)/(x * sdlog)); 30 } 31 32 33 T plnorm(T: double)(T x, T meanlog, T sdlog, int lower_tail, int log_p) 34 { 35 if (isNaN(x) || isNaN(meanlog) || isNaN(sdlog)) 36 return x + meanlog + sdlog; 37 38 if (sdlog < 0) 39 return T.nan; 40 41 if (x > 0) 42 return pnorm!T(log(x), meanlog, sdlog, lower_tail, log_p); 43 return R_DT_0!T(lower_tail, log_p); 44 } 45 46 T qlnorm(T: double)(T p, T meanlog, T sdlog, int lower_tail, int log_p) 47 { 48 if (isNaN(p) || isNaN(meanlog) || isNaN(sdlog)) 49 return p + meanlog + sdlog; 50 51 immutable T INF = T.infinity; 52 mixin (R_Q_P01_boundaries!(p, 0, INF)()); 53 54 return exp(qnorm!T(p, meanlog, sdlog, lower_tail, log_p)); 55 } 56 57 T rlnorm(T: double)(T meanlog, T sdlog) 58 { 59 if(isNaN(meanlog) || !isFinite(sdlog) || sdlog < 0.) 60 return T.nan; 61 62 return exp(rnorm!T(meanlog, sdlog)); 63 } 64 65 66 void test_lnorm() 67 { 68 import std.stdio: writeln; 69 writeln("dlnorm: ", dlnorm(10., 2., 1., 0)); 70 writeln("plnorm: ", plnorm(10., 2., 1., 1, 0)); 71 writeln("qlnorm: ", qlnorm(0.618897, 2., 1., 1, 0)); 72 writeln("rlnorm: ", rlnorm(2., 1.), ", rlnorm: ", rlnorm(2., 1.), ", rlnorm: ", rlnorm(2., 1.)); 73 } 74