R version 4.5.0 RC (2025-04-04 r88113 ucrt) -- "How About a Twenty-Six" Copyright (C) 2025 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > > ############################################################################# > ## ## > ## Tests for class 'unuran' ## > ## ## > ############################################################################# > > ## --- Test Parameters ------------------------------------------------------ > > ## size of sample for test > samplesize <- 1.e6 > > ## break points for chi^2 GoF test > nbins <- as.integer(sqrt(samplesize)) > breaks <- (0:nbins)/nbins > > ## level of significance > alpha <- 1.e-3 > > ## seed for uniform RNG > set.seed(123456) > > ## --- Load library --------------------------------------------------------- > > library(Runuran) > > > ############################################################################# > ## # > ## Auxiliary routines # > ## # > ############################################################################# > > ## -- Auxiliary routines ------------------------------------------------------ > > ## Test whether there is an error ------------------------------------------- > is.error <- function (expr) { is(try(expr), "try-error") } > > > ## --- Continuous distributions --------------------------------------------- > > ## Create an object > unr <- new("unuran", "normal()") > > ## Print object > unr Object is UNU.RAN object: method: auto distr: normal() inversion: FALSE generator ID: TDR.001 distribution: name = normal [UNU.RAN standard distribution] type = continuous univariate distribution functions = PDF dPDF domain = (-inf, inf) center = 0 [= mode] method: TDR (Transformed Density Rejection) variant = PS (proportional squeeze) T_c(x) = -1/sqrt(x) ... c = -1/2 performance characteristics: area(hat) = 1.00143 rejection constant = 1.00143 area ratio squeeze/hat = 0.994723 # intervals = 51 > print(unr) Object is UNU.RAN object: method: auto distr: normal() inversion: FALSE generator ID: TDR.001 distribution: name = normal [UNU.RAN standard distribution] type = continuous univariate distribution functions = PDF dPDF domain = (-inf, inf) center = 0 [= mode] method: TDR (Transformed Density Rejection) variant = PS (proportional squeeze) T_c(x) = -1/sqrt(x) ... c = -1/2 performance characteristics: area(hat) = 1.00143 rejection constant = 1.00143 area ratio squeeze/hat = 0.994723 # intervals = 51 > unuran.details(unr) Object is UNU.RAN object: method: auto distr: normal() inversion: FALSE generator ID: TDR.001 distribution: name = normal [UNU.RAN standard distribution] type = continuous univariate distribution functions = PDF dPDF domain = (-inf, inf) center = 0 [= mode] method: TDR (Transformed Density Rejection) variant = PS (proportional squeeze) T_c(x) = -1/sqrt(x) ... c = -1/2 performance characteristics: area(hat) = 1.00143 rejection constant = 1.00143 area ratio squeeze/hat = 0.994723 # intervals = 51 parameters: variant_ps = on [default] c = -0.5 [default] max_sqhratio = 0.99 [default] max_intervals = 100 [default] [ Hint: You may use "variant_ia" for faster generation times. ] [ Hint: You can set "max_sqhratio" closer to 1 to decrease rejection constant. ] > unuran.details(unr,show=TRUE,return.list=FALSE) Object is UNU.RAN object: method: auto distr: normal() inversion: FALSE generator ID: TDR.001 distribution: name = normal [UNU.RAN standard distribution] type = continuous univariate distribution functions = PDF dPDF domain = (-inf, inf) center = 0 [= mode] method: TDR (Transformed Density Rejection) variant = PS (proportional squeeze) T_c(x) = -1/sqrt(x) ... c = -1/2 performance characteristics: area(hat) = 1.00143 rejection constant = 1.00143 area ratio squeeze/hat = 0.994723 # intervals = 51 parameters: variant_ps = on [default] c = -0.5 [default] max_sqhratio = 0.99 [default] max_intervals = 100 [default] [ Hint: You may use "variant_ia" for faster generation times. ] [ Hint: You can set "max_sqhratio" closer to 1 to decrease rejection constant. ] > unuran.details(unr,show=TRUE,return.list=TRUE) Object is UNU.RAN object: method: auto distr: normal() inversion: FALSE generator ID: TDR.001 distribution: name = normal [UNU.RAN standard distribution] type = continuous univariate distribution functions = PDF dPDF domain = (-inf, inf) center = 0 [= mode] method: TDR (Transformed Density Rejection) variant = PS (proportional squeeze) T_c(x) = -1/sqrt(x) ... c = -1/2 performance characteristics: area(hat) = 1.00143 rejection constant = 1.00143 area ratio squeeze/hat = 0.994723 # intervals = 51 parameters: variant_ps = on [default] c = -0.5 [default] max_sqhratio = 0.99 [default] max_intervals = 100 [default] [ Hint: You may use "variant_ia" for faster generation times. ] [ Hint: You can set "max_sqhratio" closer to 1 to decrease rejection constant. ] > unuran.details(unr,show=FALSE,return.list=TRUE) > print(unuran.details(unr,show=FALSE,return.list=TRUE)) $method [1] "TDR" $type [1] "iar" $distr.class [1] "cont" $rejection.constant [1] 1.001433 $area.hat [1] 1.001433 $area.squeeze [1] 0.9961481 $intervals [1] 51 > unuran.details(unr,show=FALSE,return.list=FALSE) > print(unuran.details(unr,show=FALSE,return.list=FALSE)) NULL > > ## Test object properties > unuran.is.inversion(unr) [1] FALSE > > ## Draw samples > unuran.sample(unr) [1] 0.8342996 > unuran.sample(unr,10) [1] -0.27653690 -0.35555803 0.08741798 2.25678990 0.83502755 1.31379690 [7] 2.51009545 1.16935706 -0.42680813 -0.99728941 > x <- unuran.sample(unr, samplesize) > ur(unr) [1] 0.9941863 > ur(unr,10) [1] 0.2699989 0.3609820 -0.9613730 0.5441054 1.0018888 -0.2425041 [7] 0.5836795 -1.4984600 -0.5482153 0.3054007 > > ## Run a chi-square GoF test > chisq.test( hist(pnorm(x),plot=FALSE,breaks=breaks)$density ) Chi-squared test for given probabilities data: hist(pnorm(x), plot = FALSE, breaks = breaks)$density X-squared = 1.0234, df = 999, p-value = 1 Warning message: In chisq.test(hist(pnorm(x), plot = FALSE, breaks = breaks)$density) : Chi-squared approximation may be incorrect > > > ## Create an object > unr <- unuran.new("normal()") > > ## Draw samples > unuran.sample(unr) [1] 1.546165 > unuran.sample(unr,10) [1] 0.3659942 -1.5798172 -0.6166442 2.5150042 -1.4375299 0.7978044 [7] -0.7682332 1.1982175 1.0095702 0.3461869 > x <- unuran.sample(unr, samplesize) > > ## Run a chi-square GoF test > pval <- chisq.test( hist(pnorm(x),plot=FALSE,breaks=breaks)$density )$p.value Warning message: In chisq.test(hist(pnorm(x), plot = FALSE, breaks = breaks)$density) : Chi-squared approximation may be incorrect > if (pval < alpha) stop("chisq test FAILED! p-value=",signif(pval)) > > ## another example (for testing print) > unr <- new("unuran", "normal()", "pinv") > unr Object is UNU.RAN object: method: pinv distr: normal() inversion: TRUE generator ID: PINV.003 distribution: name = normal [UNU.RAN standard distribution] type = continuous univariate distribution functions = PDF domain = (-inf, inf) center = 0 [= mode] method: PINV (Polynomial interpolation based INVerse CDF) order of polynomial = 5 smoothness = 0 [continuous] use PDF + Lobatto integration [default] performance characteristics: truncated domain = (-6.82147,6.82147) u-error <= 8.78461e-11 (mean = 2.92896e-11) [ u-resolution = 1e-10 ] area below PDF = 0.99999999999102174 # intervals = 124 > print(unr) Object is UNU.RAN object: method: pinv distr: normal() inversion: TRUE generator ID: PINV.003 distribution: name = normal [UNU.RAN standard distribution] type = continuous univariate distribution functions = PDF domain = (-inf, inf) center = 0 [= mode] method: PINV (Polynomial interpolation based INVerse CDF) order of polynomial = 5 smoothness = 0 [continuous] use PDF + Lobatto integration [default] performance characteristics: truncated domain = (-6.82147,6.82147) u-error <= 8.78461e-11 (mean = 2.92896e-11) [ u-resolution = 1e-10 ] area below PDF = 0.99999999999102174 # intervals = 124 > unuran.details(unr) Object is UNU.RAN object: method: pinv distr: normal() inversion: TRUE generator ID: PINV.003 distribution: name = normal [UNU.RAN standard distribution] type = continuous univariate distribution functions = PDF domain = (-inf, inf) center = 0 [= mode] method: PINV (Polynomial interpolation based INVerse CDF) order of polynomial = 5 smoothness = 0 [continuous] use PDF + Lobatto integration [default] performance characteristics: truncated domain = (-6.82147,6.82147) u-error <= 8.78461e-11 (mean = 2.92896e-11) [ u-resolution = 1e-10 ] area below PDF = 0.99999999999102174 # intervals = 124 parameters: order = 5 [default] smoothness = 0 [default] u_resolution = 1e-10 [default] use_upoints = FALSE [default] boundary = (-1e+100,1e+100) [default] search for boundary: left=TRUE, right=TRUE [default] maximum number of interval = 10000 [default] keep table of CDF values = FALSE [default] [ Hint: You can increase "order" to decrease #intervals ] [ Hint: You can decrease the u-error by decreasing "u_resolution". (it is bounded by the machine epsilon, however.) ] > print(unuran.details(unr,show=FALSE,return.list=TRUE)) $method [1] "PINV" $type [1] "inv" $distr.class [1] "cont" $truncated.domain [1] -6.82147 6.82147 $area.pdf [1] 1 $intervals [1] 124 > print(unuran.details(unr,show=FALSE,debug=TRUE)) $method [1] "PINV" $type [1] "inv" $distr.class [1] "cont" $truncated.domain [1] -6.82147 6.82147 $area.pdf [1] 1 $intervals [1] 124 $cdfi [1] 0.000000e+00 4.905094e-12 3.519950e-11 3.676356e-10 2.519630e-09 [6] 1.110170e-08 3.996927e-08 1.212936e-07 3.182410e-07 7.373627e-07 [11] 1.535796e-06 3.131404e-06 5.838675e-06 1.070114e-05 1.819156e-05 [16] 3.049949e-05 4.798875e-05 7.466777e-05 1.100999e-04 1.608892e-04 [21] 2.330031e-04 3.344261e-04 4.594327e-04 6.266184e-04 8.485000e-04 [26] 1.140713e-03 1.479736e-03 1.908426e-03 2.447128e-03 3.119863e-03 [31] 3.954755e-03 4.984447e-03 6.246494e-03 7.783709e-03 9.644442e-03 [36] 1.188278e-02 1.426950e-02 1.705814e-02 2.029995e-02 2.404959e-02 [41] 2.836474e-02 3.330566e-02 3.893458e-02 4.531499e-02 5.251074e-02 [46] 6.058512e-02 6.959976e-02 7.961341e-02 9.068072e-02 1.028509e-01 [51] 1.161664e-01 1.337045e-01 1.529821e-01 1.740182e-01 1.968067e-01 [56] 2.213148e-01 2.474809e-01 2.809399e-01 3.164443e-01 3.537267e-01 [61] 3.924682e-01 4.403713e-01 4.891644e-01 5.478664e-01 6.039486e-01 [66] 6.526773e-01 6.944972e-01 7.300984e-01 7.634926e-01 7.914914e-01 [71] 8.174346e-01 8.412616e-01 8.608794e-01 8.787782e-01 8.949923e-01 [76] 9.095755e-01 9.225982e-01 9.330538e-01 9.423819e-01 9.506560e-01 [81] 9.579528e-01 9.643504e-01 9.699274e-01 9.747608e-01 9.789256e-01 [86] 9.824935e-01 9.855324e-01 9.881057e-01 9.902723e-01 9.920858e-01 [91] 9.935950e-01 9.948437e-01 9.958709e-01 9.967110e-01 9.973941e-01 [96] 9.979464e-01 9.983903e-01 9.987450e-01 9.990268e-01 9.992494e-01 [101] 9.994242e-01 9.995607e-01 9.996848e-01 9.997756e-01 9.998415e-01 [106] 9.998889e-01 9.999228e-01 9.999467e-01 9.999662e-01 9.999788e-01 [111] 9.999869e-01 9.999920e-01 9.999956e-01 9.999976e-01 9.999988e-01 [116] 9.999994e-01 9.999997e-01 9.999999e-01 1.000000e+00 1.000000e+00 [121] 1.000000e+00 1.000000e+00 1.000000e+00 1.000000e+00 1.000000e+00 $xi [1] -6.82146959 -6.71488413 -6.50171321 -6.15637631 -5.84557310 -5.59382250 [7] -5.36724696 -5.16332898 -4.97980279 -4.81462922 -4.66597301 -4.51731680 [13] -4.38352621 -4.24973562 -4.12932409 -4.00891256 -3.90054218 -3.79217180 [19] -3.69463846 -3.59710512 -3.49957178 -3.40203844 -3.31425844 -3.22647843 [25] -3.13869842 -3.05091842 -2.97191641 -2.89291441 -2.81391240 -2.73491040 [31] -2.65590839 -2.57690639 -2.49790438 -2.41890237 -2.33990037 -2.26089836 [37] -2.18979656 -2.11869475 -2.04759295 -1.97649114 -1.90538934 -1.83428753 [43] -1.76318573 -1.69208392 -1.62098212 -1.54988031 -1.47877851 -1.40767670 [49] -1.33657490 -1.26547309 -1.19437129 -1.10904912 -1.02372696 -0.93840479 [55] -0.85308263 -0.76776046 -0.68243829 -0.58005170 -0.47766510 -0.37527850 [61] -0.27289190 -0.15002798 -0.02716406 0.12027264 0.26358112 0.39255874 [67] 0.50863861 0.61311049 0.71758237 0.81160706 0.90563175 0.99965644 [73] 1.08427866 1.16890088 1.25352310 1.33814532 1.42276754 1.49892754 [79] 1.57508754 1.65124754 1.72740754 1.80356754 1.87972754 1.95588754 [85] 2.03204754 2.10820754 2.18436754 2.26052754 2.33668754 2.41284753 [91] 2.48900753 2.56516753 2.64132753 2.71748753 2.79364753 2.86980753 [97] 2.94596753 3.02212753 3.09828753 3.17444753 3.25060753 3.32676753 [103] 3.41815952 3.50955152 3.60094352 3.69233552 3.78372752 3.87511952 [109] 3.98478992 4.09446032 4.20413072 4.31380112 4.44540559 4.57701007 [115] 4.70861455 4.86653993 5.02446530 5.18239068 5.37190112 5.56141157 [121] 5.78882411 6.06171916 6.33461421 6.65404003 6.82146959 $ui [,1] [,2] [,3] [,4] [,5] [1,] 4.704858e-13 1.401776e-12 2.690599e-12 4.026347e-12 4.905094e-12 [2,] 2.031467e-12 6.611201e-12 1.415293e-11 2.341424e-11 3.029441e-11 [3,] 1.424484e-11 5.159075e-11 1.265823e-10 2.379088e-10 3.324361e-10 [4,] 1.114417e-10 3.860423e-10 8.941911e-10 1.591673e-09 2.151995e-09 [5,] 5.652676e-10 1.848501e-09 3.979490e-09 6.615187e-09 8.582066e-09 [6,] 2.117132e-09 6.744978e-09 1.405332e-08 2.266248e-08 2.886758e-08 [7,] 6.518197e-09 2.031840e-08 4.120326e-08 6.480750e-08 8.132434e-08 [8,] 1.699444e-08 5.201373e-08 1.031283e-07 1.588970e-07 1.969474e-07 [9,] 3.846571e-08 1.159259e-07 2.255572e-07 3.416233e-07 4.191218e-07 [10,] 7.717193e-08 2.295561e-07 4.396263e-07 6.563682e-07 7.984329e-07 [11,] 1.558989e-07 4.624810e-07 8.827841e-07 1.314057e-06 1.595608e-06 [12,] 2.759918e-07 8.099271e-07 1.525992e-06 2.245065e-06 2.707271e-06 [13,] 4.999863e-07 1.464069e-06 2.751215e-06 4.038018e-06 4.862469e-06 [14,] 7.976629e-07 2.314807e-06 4.303369e-06 6.255744e-06 7.490414e-06 [15,] 1.319716e-06 3.823091e-06 7.092427e-06 1.029069e-05 1.230794e-05 [16,] 1.930510e-06 5.550673e-06 1.020591e-05 1.469117e-05 1.748925e-05 [17,] 2.961114e-06 8.501896e-06 1.560596e-05 2.243061e-05 2.667902e-05 [18,] 4.028814e-06 1.149496e-05 2.094371e-05 2.990509e-05 3.543217e-05 [19,] 5.800456e-06 1.653097e-05 3.007861e-05 4.289701e-05 5.078929e-05 [20,] 8.272099e-06 2.354825e-05 4.278899e-05 6.095089e-05 7.211383e-05 [21,] 1.168525e-05 3.322675e-05 6.029429e-05 8.578358e-05 1.014231e-04 [22,] 1.468256e-05 4.153928e-05 7.493237e-05 1.060539e-04 1.250065e-04 [23,] 1.970592e-05 5.570012e-05 1.003683e-04 1.419180e-04 1.671857e-04 [24,] 2.624493e-05 7.411519e-05 1.334068e-04 1.884531e-04 2.218816e-04 [25,] 3.468548e-05 9.786159e-05 1.759601e-04 2.483272e-04 2.922128e-04 [26,] 4.086622e-05 1.148330e-04 2.054975e-04 2.888073e-04 3.390234e-04 [27,] 5.182079e-05 1.455076e-04 2.601645e-04 3.653565e-04 4.286903e-04 [28,] 6.530300e-05 1.832290e-04 3.273252e-04 4.593206e-04 5.387017e-04 [29,] 8.178085e-05 2.292938e-04 4.092616e-04 5.738594e-04 6.727348e-04 [30,] 1.017793e-04 2.851544e-04 5.085252e-04 7.125011e-04 8.348919e-04 [31,] 1.258801e-04 3.524175e-04 6.279340e-04 8.791361e-04 1.029692e-03 [32,] 1.547191e-04 4.328371e-04 7.705583e-04 1.077996e-03 1.262047e-03 [33,] 1.889819e-04 5.283006e-04 9.396951e-04 1.313617e-03 1.537214e-03 [34,] 2.293961e-04 6.408071e-04 1.138829e-03 1.590782e-03 1.860734e-03 [35,] 2.767205e-04 7.724372e-04 1.371576e-03 1.914446e-03 2.238333e-03 [36,] 2.982006e-04 8.300583e-04 1.469085e-03 2.044727e-03 2.386721e-03 [37,] 3.492027e-04 9.714483e-04 1.718130e-03 2.389902e-03 2.788639e-03 [38,] 4.068659e-04 1.131190e-03 1.999264e-03 2.779265e-03 3.241814e-03 [39,] 4.716603e-04 1.310556e-03 2.314670e-03 3.215769e-03 3.749638e-03 [40,] 5.440163e-04 1.510708e-03 2.666324e-03 3.702072e-03 4.315151e-03 [41,] 6.243080e-04 1.732646e-03 3.055916e-03 4.240432e-03 4.940922e-03 [42,] 7.128372e-04 1.977168e-03 3.484775e-03 4.832595e-03 5.628923e-03 [43,] 8.098159e-04 2.244823e-03 3.953784e-03 5.479688e-03 6.380402e-03 [44,] 9.153489e-04 2.535859e-03 4.463298e-03 6.182106e-03 7.195750e-03 [45,] 1.029417e-03 2.850183e-03 5.013069e-03 6.939404e-03 8.074386e-03 [46,] 1.151863e-03 3.187315e-03 5.602169e-03 7.750202e-03 9.014638e-03 [47,] 1.282373e-03 3.546351e-03 6.228933e-03 8.612099e-03 1.001365e-02 [48,] 1.420472e-03 3.925934e-03 6.890898e-03 9.521604e-03 1.106731e-02 [49,] 1.565508e-03 4.324230e-03 7.584778e-03 1.047409e-02 1.217018e-02 [50,] 1.716653e-03 4.738919e-03 8.306435e-03 1.146378e-02 1.331549e-02 [51,] 2.250026e-03 6.219726e-03 1.091862e-02 1.508823e-02 1.753811e-02 [52,] 2.481143e-03 6.852764e-03 1.201799e-02 1.659308e-02 1.927763e-02 [53,] 2.716155e-03 7.495472e-03 1.313213e-02 1.811572e-02 2.103608e-02 [54,] 2.951859e-03 8.138995e-03 1.424550e-02 1.963469e-02 2.278853e-02 [55,] 3.184749e-03 8.773669e-03 1.534121e-02 2.112673e-02 2.450801e-02 [56,] 3.411091e-03 9.389238e-03 1.640138e-02 2.256734e-02 2.616617e-02 [57,] 4.355781e-03 1.199492e-02 2.096240e-02 2.885242e-02 3.345895e-02 [58,] 4.643330e-03 1.277114e-02 2.228730e-02 3.063816e-02 3.550442e-02 [59,] 4.898244e-03 1.345582e-02 2.344892e-02 3.219537e-02 3.728238e-02 [60,] 5.113269e-03 1.402937e-02 2.441389e-02 3.347914e-02 3.874148e-02 [61,] 6.340771e-03 1.738583e-02 3.022859e-02 4.141806e-02 4.790311e-02 [62,] 6.501087e-03 1.779415e-02 3.087582e-02 4.223051e-02 4.879317e-02 [63,] 7.878911e-03 2.152654e-02 3.726888e-02 5.086976e-02 5.870197e-02 [64,] 7.595131e-03 2.070065e-02 3.573934e-02 4.866587e-02 5.608223e-02 [65,] 6.642788e-03 1.807130e-02 3.113560e-02 4.232487e-02 4.872865e-02 [66,] 5.726431e-03 1.555852e-02 2.676879e-02 3.634710e-02 4.181993e-02 [67,] 4.888706e-03 1.327149e-02 2.281355e-02 3.095432e-02 3.560114e-02 [68,] 4.607115e-03 1.249130e-02 2.144137e-02 2.905626e-02 3.339423e-02 [69,] 3.867116e-03 1.048122e-02 1.798456e-02 2.436515e-02 2.799885e-02 [70,] 3.596736e-03 9.738469e-03 1.669055e-02 2.258931e-02 2.594317e-02 [71,] 3.315818e-03 8.968717e-03 1.535334e-02 2.075867e-02 2.382700e-02 [72,] 2.728669e-03 7.381294e-03 1.263773e-02 1.708970e-02 1.961777e-02 [73,] 2.497179e-03 6.749524e-03 1.154515e-02 1.559956e-02 1.789885e-02 [74,] 2.269021e-03 6.127793e-03 1.047178e-02 1.413779e-02 1.621409e-02 [75,] 2.046998e-03 5.523640e-03 9.430452e-03 1.272161e-02 1.458316e-02 [76,] 1.833524e-03 4.943526e-03 8.432089e-03 1.136565e-02 1.302276e-02 [77,] 1.468721e-03 3.962274e-03 6.763121e-03 9.121750e-03 1.045556e-02 [78,] 1.313567e-03 3.541341e-03 6.040034e-03 8.141164e-03 9.328095e-03 [79,] 1.168008e-03 3.146822e-03 5.363065e-03 7.223982e-03 8.274100e-03 [80,] 1.032573e-03 2.780083e-03 4.734435e-03 6.373070e-03 7.296771e-03 [81,] 9.075626e-04 2.441881e-03 4.155322e-03 5.589881e-03 6.397684e-03 [82,] 7.930733e-04 2.132417e-03 3.625957e-03 4.874593e-03 5.576954e-03 [83,] 6.890187e-04 1.851403e-03 3.145734e-03 4.226257e-03 4.833408e-03 [84,] 5.951544e-04 1.598126e-03 2.713331e-03 3.642968e-03 4.164780e-03 [85,] 5.111040e-04 1.371520e-03 2.326832e-03 3.122027e-03 3.567901e-03 [86,] 4.363851e-04 1.170238e-03 1.983850e-03 2.660112e-03 3.038895e-03 [87,] 3.704345e-04 9.927218e-04 1.681643e-03 2.253434e-03 2.573361e-03 [88,] 3.126324e-04 8.372631e-04 1.417230e-03 1.897893e-03 2.166546e-03 [89,] 2.623236e-04 7.020652e-04 1.187486e-03 1.589207e-03 1.813499e-03 [90,] 2.188376e-04 5.852940e-04 9.892311e-04 1.323034e-03 1.509207e-03 [91,] 1.815044e-04 4.851229e-04 8.193107e-04 1.095074e-03 1.248713e-03 [92,] 1.496696e-04 3.997704e-04 6.746536e-04 9.011513e-04 1.027208e-03 [93,] 1.227046e-04 3.275296e-04 5.523248e-04 7.372824e-04 8.401104e-04 [94,] 1.000159e-04 2.667912e-04 4.495621e-04 5.997246e-04 6.831192e-04 [95,] 8.105103e-05 2.160596e-04 3.638030e-04 4.850111e-04 5.522539e-04 [96,] 6.530234e-05 1.739630e-04 2.927010e-04 3.899720e-04 4.438776e-04 [97,] 5.230942e-05 1.392583e-04 2.341336e-04 3.117432e-04 3.547071e-04 [98,] 4.165931e-05 1.108324e-04 1.862021e-04 2.477664e-04 2.818115e-04 [99,] 3.298566e-05 8.769871e-05 1.472268e-04 1.957806e-04 2.226023e-04 [100,] 2.596685e-05 6.899234e-05 1.157365e-04 1.538079e-04 1.748166e-04 [101,] 2.032330e-05 5.396219e-05 9.045564e-05 1.201350e-04 1.364954e-04 [102,] 1.891312e-05 4.989509e-05 8.302779e-05 1.095850e-04 1.240675e-04 [103,] 1.388882e-05 3.660580e-05 6.084919e-05 8.024067e-05 9.079917e-05 [104,] 1.011440e-05 2.663268e-05 4.422417e-05 5.826561e-05 6.589932e-05 [105,] 7.304448e-06 1.921555e-05 3.187411e-05 4.195700e-05 4.743020e-05 [106,] 5.231272e-06 1.374876e-05 2.278191e-05 2.996202e-05 3.385354e-05 [107,] 3.715350e-06 9.755447e-06 1.614789e-05 2.121837e-05 2.396229e-05 [108,] 3.125351e-06 8.133083e-06 1.332751e-05 1.736448e-05 1.951626e-05 [109,] 2.029418e-06 5.274044e-06 8.629554e-06 1.122952e-05 1.261228e-05 [110,] 1.302029e-06 3.379165e-06 5.520855e-06 7.175315e-06 8.053271e-06 [111,] 8.253658e-07 2.139204e-06 3.489819e-06 4.530030e-06 5.080807e-06 [112,] 6.164569e-07 1.578789e-06 2.541635e-06 3.263035e-06 3.637256e-06 [113,] 3.460100e-07 8.844624e-07 1.420883e-06 1.821047e-06 2.027966e-06 [114,] 1.908767e-07 4.869836e-07 7.807035e-07 9.988653e-07 1.111313e-06 [115,] 1.231733e-07 3.093868e-07 4.875898e-07 6.151892e-07 6.791860e-07 [116,] 5.773537e-08 1.446272e-07 2.272667e-07 2.860700e-07 3.154264e-07 [117,] 2.639589e-08 6.594315e-08 1.033224e-07 1.297540e-07 1.428886e-07 [118,] 1.397336e-08 3.420204e-08 5.242726e-08 6.469710e-08 7.057447e-08 [119,] 5.128054e-09 1.250390e-08 1.909021e-08 2.348448e-08 2.557535e-08 [120,] 2.148801e-09 5.104198e-09 7.582866e-09 9.133224e-09 9.836020e-09 [121,] 6.949912e-10 1.597942e-09 2.296921e-09 2.699322e-09 2.870672e-09 [122,] 1.372925e-10 3.133162e-10 4.470888e-10 5.226676e-10 5.544032e-10 [123,] 2.891352e-11 6.360547e-11 8.761100e-11 9.991289e-11 1.047102e-10 [124,] 2.019558e-12 4.875282e-12 7.367024e-12 8.991762e-12 9.752061e-12 $zi [,1] [,2] [,3] [,4] [,5] [1,] 3.035106e+10 -2.705859e+21 2.692938e+32 -2.544437e+43 2.301958e+54 [2,] 1.405855e+10 -4.927063e+20 1.600006e+31 -4.169121e+41 9.391915e+51 [3,] 3.247940e+09 -2.136373e+19 1.046480e+29 -3.391740e+38 8.490063e+47 [4,] 3.736459e+08 -2.875389e+17 1.804689e+26 -8.112697e+34 2.949378e+43 [5,] 5.966765e+07 -7.650995e+15 9.017820e+23 -8.456342e+31 6.822133e+39 [6,] 1.433796e+07 -4.413574e+14 1.319967e+22 -3.297671e+29 7.296940e+36 [7,] 4.191317e+06 -3.752410e+13 3.420736e+20 -2.714809e+27 1.956307e+34 [8,] 1.446817e+06 -4.435449e+12 1.440944e+19 -4.221024e+25 1.146845e+32 [9,] 5.752933e+05 -6.941976e+11 9.184899e+17 -1.129111e+24 1.311134e+30 [10,] 2.580751e+05 -1.380978e+11 8.343220e+16 -4.804943e+22 2.655344e+28 [11,] 1.277505e+05 -3.293798e+10 9.654732e+15 -2.712393e+21 7.336367e+26 [12,] 6.494591e+04 -8.395688e+09 1.264802e+15 -1.865977e+20 2.685779e+25 [13,] 3.585006e+04 -2.491140e+09 2.029293e+14 -1.626095e+19 1.274696e+24 [14,] 2.022419e+04 -7.806076e+08 3.607389e+13 -1.669882e+18 7.647472e+22 [15,] 1.222391e+04 -2.778751e+08 7.602072e+12 -2.090950e+17 5.702476e+21 [16,] 7.520745e+03 -1.034563e+08 1.744030e+12 -3.001555e+16 5.171192e+20 [17,] 4.903181e+03 -4.288148e+07 4.616876e+11 -5.090311e+15 5.628483e+19 [18,] 3.243384e+03 -1.844354e+07 1.311498e+11 -9.675455e+14 7.216983e+18 [19,] 2.252752e+03 -8.684512e+06 4.203761e+10 -2.116470e+14 1.079015e+18 [20,] 1.579646e+03 -4.164888e+06 1.384685e+10 -4.800529e+13 1.687843e+17 [21,] 1.118246e+03 -2.034235e+06 4.686832e+09 -1.128926e+13 2.761994e+16 [22,] 8.009699e+02 -1.023658e+06 1.679418e+09 -2.911598e+12 5.161585e+15 [23,] 5.967896e+02 -5.544186e+05 6.638762e+08 -8.418661e+11 1.093030e+15 [24,] 4.480977e+02 -3.047234e+05 2.682344e+08 -2.505968e+11 2.400085e+14 [25,] 3.390552e+02 -1.699580e+05 1.107683e+08 -7.678822e+10 5.464168e+13 [26,] 2.589978e+02 -9.710092e+04 4.789926e+07 -2.535672e+10 1.385510e+13 [27,] 2.042474e+02 -5.889059e+04 2.243614e+07 -9.189979e+09 3.889566e+12 [28,] 1.620793e+02 -3.613929e+04 1.069482e+07 -3.409120e+09 1.124093e+12 [29,] 1.294222e+02 -2.243920e+04 5.187764e+06 -1.294330e+09 3.344057e+11 [30,] 1.039922e+02 -1.409647e+04 2.560609e+06 -5.029092e+08 1.023941e+11 [31,] 8.408210e+01 -8.959177e+03 1.285982e+06 -1.999588e+08 3.226728e+10 [32,] 6.840954e+01 -5.760479e+03 6.570929e+05 -8.135059e+07 1.046383e+10 [33,] 5.600674e+01 -3.746774e+03 3.415763e+05 -3.386186e+07 3.491497e+09 [34,] 4.613968e+01 -2.465120e+03 1.806280e+05 -1.441946e+07 1.198600e+09 [35,] 3.824893e+01 -1.640480e+03 9.715925e+04 -6.281038e+06 4.232771e+08 [36,] 3.194439e+01 -1.111373e+03 5.406254e+04 -2.889339e+06 1.616325e+08 [37,] 2.727881e+01 -7.856165e+02 3.185880e+04 -1.421578e+06 6.645789e+07 [38,] 2.341272e+01 -5.603884e+02 1.903124e+04 -7.122015e+05 2.794988e+07 [39,] 2.019639e+01 -4.033330e+02 1.152334e+04 -3.632912e+05 1.202211e+07 [40,] 1.751020e+01 -2.928865e+02 7.071794e+03 -1.886622e+05 5.288087e+06 [41,] 1.525823e+01 -2.145639e+02 4.398322e+03 -9.973527e+04 2.378374e+06 [42,] 1.336327e+01 -1.585595e+02 2.772136e+03 -5.366613e+04 1.093629e+06 [43,] 1.176296e+01 -1.181836e+02 1.770417e+03 -2.938940e+04 5.140549e+05 [44,] 1.040678e+01 -8.883758e+01 1.145599e+03 -1.637830e+04 2.469663e+05 [45,] 9.253619e+00 -6.733624e+01 7.510139e+02 -9.287114e+03 1.212524e+05 [46,] 8.269939e+00 -5.145686e+01 4.987512e+02 -5.357591e+03 6.082755e+04 [47,] 7.428286e+00 -3.963689e+01 3.355094e+02 -3.143946e+03 3.117437e+04 [48,] 6.706106e+00 -3.076984e+01 2.286004e+02 -1.876434e+03 1.631959e+04 [49,] 6.084820e+00 -2.406654e+01 1.577510e+02 -1.138873e+03 8.724868e+03 [50,] 5.549076e+00 -1.896008e+01 1.102475e+02 -7.027912e+02 4.762862e+03 [51,] 5.080388e+00 -1.492643e+01 7.654651e+01 -4.250754e+02 2.507134e+03 [52,] 4.607152e+00 -1.140430e+01 5.145442e+01 -2.485346e+02 1.280233e+03 [53,] 4.208524e+00 -8.787369e+00 3.522796e+01 -1.486706e+02 6.730735e+02 [54,] 3.872475e+00 -6.821165e+00 2.457197e+01 -9.093474e+01 3.642287e+02 [55,] 3.589294e+00 -5.326766e+00 1.747020e+01 -5.682788e+01 2.028321e+02 [56,] 3.351128e+00 -4.177105e+00 1.267114e+01 -3.624374e+01 1.162358e+02 [57,] 3.149195e+00 -3.254731e+00 9.236069e+00 -2.277586e+01 6.518486e+01 [58,] 2.954174e+00 -2.428460e+00 6.655225e+00 -1.387953e+01 3.614371e+01 [59,] 2.800433e+00 -1.789142e+00 4.981132e+00 -8.576100e+00 2.103642e+01 [60,] 2.682668e+00 -1.279546e+00 3.897100e+00 -5.280585e+00 1.299400e+01 [61,] 2.596000e+00 -8.452567e-01 3.183671e+00 -2.974360e+00 8.418468e+00 [62,] 2.531983e+00 -4.139371e-01 2.750456e+00 -1.197683e+00 6.208085e+00 [63,] 2.507044e+00 -8.142855e-03 2.629225e+00 3.681544e-01 5.908910e+00 [64,] 2.527895e+00 4.624735e-01 2.872537e+00 2.241612e+00 7.919183e+00 [65,] 2.601276e+00 9.715512e-01 3.548406e+00 4.811457e+00 1.323639e+01 [66,] 2.715785e+00 1.533475e+00 4.693034e+00 8.670570e+00 2.379555e+01 [67,] 2.863044e+00 2.179678e+00 6.405083e+00 1.464899e+01 4.325395e+01 [68,] 3.038035e+00 2.948069e+00 8.911157e+00 2.434121e+01 8.011955e+01 [69,] 3.257446e+00 3.943916e+00 1.270250e+01 4.053669e+01 1.499073e+02 [70,] 3.502319e+00 5.152823e+00 1.812034e+01 6.698699e+01 2.814034e+02 [71,] 3.799038e+00 6.763862e+00 2.647930e+01 1.135271e+02 5.469357e+02 [72,] 4.154856e+00 8.900656e+00 3.913465e+01 1.931235e+02 1.061880e+03 [73,] 4.540014e+00 1.153083e+01 5.727379e+01 3.252074e+02 2.050739e+03 [74,] 4.996528e+00 1.506375e+01 8.534299e+01 5.602013e+02 4.076568e+03 [75,] 5.538464e+00 1.986083e+01 1.294672e+02 9.875499e+02 8.343843e+03 [76,] 6.183300e+00 2.644488e+01 1.999556e+02 1.782180e+03 1.758980e+04 [77,] 6.947205e+00 3.540117e+01 3.103381e+02 3.211903e+03 3.672028e+04 [78,] 7.767786e+00 4.665984e+01 4.734546e+02 5.692812e+03 7.562567e+04 [79,] 8.735814e+00 6.205873e+01 7.330339e+02 1.028984e+04 1.596298e+05 [80,] 9.881630e+00 8.331058e+01 1.151904e+03 1.897059e+04 3.454021e+05 [81,] 1.124276e+01 1.129080e+02 1.837387e+03 3.567901e+04 7.662726e+05 [82,] 1.286578e+01 1.545092e+02 2.975253e+03 6.846516e+04 1.743274e+06 [83,] 1.480875e+01 2.135293e+02 4.891387e+03 1.340643e+05 4.067670e+06 [84,] 1.714430e+01 2.980539e+02 8.165256e+03 2.679182e+05 9.736271e+06 [85,] 1.996366e+01 4.202614e+02 1.384143e+04 5.465024e+05 2.390963e+07 [86,] 2.338188e+01 5.986586e+02 2.382905e+04 1.137982e+06 6.024888e+07 [87,] 2.754470e+01 8.616216e+02 4.166653e+04 2.419257e+06 1.558048e+08 [88,] 3.263739e+01 1.253054e+03 7.400493e+04 5.251450e+06 4.135469e+08 [89,] 3.889663e+01 1.841505e+03 1.335249e+05 1.164051e+07 1.126769e+09 [90,] 4.662593e+01 2.735001e+03 2.447522e+05 2.635132e+07 3.151821e+09 [91,] 5.621628e+01 4.105364e+03 4.558116e+05 6.092713e+07 9.052188e+09 [92,] 6.817353e+01 6.228466e+03 8.625195e+05 1.438915e+08 2.669664e+10 [93,] 8.315501e+01 9.551449e+03 1.658463e+06 3.471454e+08 8.085640e+10 [94,] 1.020188e+02 1.480599e+04 3.240582e+06 8.556057e+08 2.515174e+11 [95,] 1.258899e+02 2.320095e+04 6.434961e+06 2.154538e+09 8.036311e+11 [96,] 1.562502e+02 3.675292e+04 1.298669e+07 5.543495e+09 2.637654e+12 [97,] 1.950606e+02 5.885898e+04 2.663808e+07 1.457440e+10 8.893782e+12 [98,] 2.449274e+02 9.529796e+04 5.553679e+07 3.915644e+10 3.081034e+13 [99,] 3.093315e+02 1.559979e+05 1.176934e+08 1.075095e+11 1.096679e+14 [100,] 3.929435e+02 2.581854e+05 2.535338e+08 3.016799e+11 4.011123e+14 [101,] 5.020594e+02 4.320504e+05 5.552004e+08 8.652171e+11 1.507595e+15 [102,] 6.473920e+02 7.440006e+05 1.293475e+09 2.764901e+12 6.669944e+15 [103,] 8.815872e+02 1.419710e+06 3.457515e+09 1.036668e+13 3.511260e+16 [104,] 1.210572e+03 2.752788e+06 9.464809e+09 4.011882e+13 1.922935e+17 [105,] 1.676267e+03 5.423859e+06 2.653539e+10 1.602645e+14 1.095638e+18 [106,] 2.340579e+03 1.085978e+07 7.619526e+10 6.609034e+14 6.495433e+18 [107,] 3.295572e+03 2.209654e+07 2.240993e+11 2.813704e+15 4.007039e+19 [108,] 4.701248e+03 4.681049e+07 7.186024e+11 1.392102e+16 3.098813e+20 [109,] 7.240032e+03 1.144153e+08 2.790624e+12 8.608116e+16 3.055855e+21 [110,] 1.128473e+04 2.862500e+08 1.121997e+13 5.574335e+17 3.192232e+22 [111,] 1.780186e+04 7.330678e+08 4.670766e+13 3.780631e+18 3.532880e+23 [112,] 2.860160e+04 1.984242e+09 2.186131e+14 3.139181e+19 5.295924e+24 [113,] 5.095707e+04 6.511357e+09 1.324662e+15 3.524233e+20 1.104172e+26 [114,] 9.237196e+04 2.210108e+10 8.441115e+15 4.230512e+21 2.502844e+27 [115,] 1.717741e+05 8.088697e+10 6.293857e+16 6.646636e+22 8.480944e+28 [116,] 3.664649e+05 3.822664e+11 6.617839e+17 1.562801e+24 4.474910e+30 [117,] 8.015638e+05 1.896972e+12 7.484383e+18 4.048503e+25 2.664847e+32 [118,] 1.817000e+06 1.043641e+13 1.046297e+20 1.503777e+27 2.713673e+34 [119,] 4.951115e+06 8.085908e+13 2.321103e+21 9.622985e+28 5.035685e+36 [120,] 1.417884e+07 7.203978e+14 6.827183e+22 9.901511e+30 1.889408e+39 [121,] 5.260642e+07 1.096742e+16 4.596243e+24 3.172331e+33 3.040721e+42 [122,] 2.663000e+08 2.982792e+17 6.819369e+26 2.607857e+36 1.401074e+46 [123,] 1.480101e+09 1.032152e+19 1.600316e+29 4.517808e+39 1.911627e+50 [124,] 1.110704e+10 5.123757e+20 4.212696e+31 5.085974e+42 7.841903e+53 > unuran.is.inversion(unr) [1] TRUE > > ## remove (so that valgrind does not see lost memory from UNU.RAN) > rm(unr) > > > ## --- Continuous distributions - S4 distribution object -------------------- > > ## use PDF > gausspdf <- function (x) { exp(-0.5*x^2) } > gaussdpdf <- function (x) { -x*exp(-0.5*x^2) } > gauss <- new("unuran.cont", pdf=gausspdf, dpdf=gaussdpdf, lb=-Inf, ub=Inf, center=0.1) > unr <- 0; unr <- unuran.new(gauss, "tdr") > unr Object is UNU.RAN object: method: tdr distr: [S4 class] inversion: FALSE generator ID: TDR.004 distribution: name = unknown type = continuous univariate distribution functions = PDF dPDF domain = (-inf, inf) center = 0.1 method: TDR (Transformed Density Rejection) variant = PS (proportional squeeze) T_c(x) = -1/sqrt(x) ... c = -1/2 performance characteristics: area(hat) = 2.50984 rejection constant <= 1.00467 area ratio squeeze/hat = 0.995355 # intervals = 54 > x <- unuran.sample(unr, samplesize) > pval <- chisq.test( hist(pnorm(x),plot=FALSE,breaks=breaks)$density )$p.value Warning message: In chisq.test(hist(pnorm(x), plot = FALSE, breaks = breaks)$density) : Chi-squared approximation may be incorrect > if (pval < alpha) stop("chisq test FAILED! p-value=",signif(pval)) > rm(unr) > > ## use PDF > gausspdf <- function (x) { exp(-0.5*x^2) } > gaussdpdf <- function (x) { -x*exp(-0.5*x^2) } > gauss <- unuran.cont.new(pdf=gausspdf, dpdf=gaussdpdf, lb=-Inf, ub=Inf, center=0.1) > unr <- unuran.new(gauss, "tdr") > unr Object is UNU.RAN object: method: tdr distr: [S4 class] inversion: FALSE generator ID: TDR.005 distribution: name = unknown type = continuous univariate distribution functions = PDF dPDF domain = (-inf, inf) center = 0.1 method: TDR (Transformed Density Rejection) variant = PS (proportional squeeze) T_c(x) = -1/sqrt(x) ... c = -1/2 performance characteristics: area(hat) = 2.50984 rejection constant <= 1.00467 area ratio squeeze/hat = 0.995355 # intervals = 54 > x <- unuran.sample(unr, samplesize) > pval <- chisq.test( hist(pnorm(x),plot=FALSE,breaks=breaks)$density )$p.value Warning message: In chisq.test(hist(pnorm(x), plot = FALSE, breaks = breaks)$density) : Chi-squared approximation may be incorrect > if (pval < alpha) stop("chisq test FAILED! p-value=",signif(pval)) > rm(unr) > > ## use logPDF > gausspdf <- function (x) { -0.5*x^2 } > gaussdpdf <- function (x) { -x } > gauss <- new("unuran.cont", pdf=gausspdf, dpdf=gaussdpdf, islog=TRUE, lb=-Inf, ub=Inf, mode=0) > unr <- unuran.new(gauss, "tdr") > unr Object is UNU.RAN object: method: tdr distr: [S4 class] inversion: FALSE generator ID: TDR.006 distribution: name = unknown type = continuous univariate distribution functions = PDF dPDF domain = (-inf, inf) center = 0 [= mode] method: TDR (Transformed Density Rejection) variant = PS (proportional squeeze) T_c(x) = -1/sqrt(x) ... c = -1/2 performance characteristics: area(hat) = 2.51022 rejection constant <= 1.00531 area ratio squeeze/hat = 0.994723 # intervals = 51 > x <- unuran.sample(unr, samplesize) > pval <- chisq.test( hist(pnorm(x),plot=FALSE,breaks=breaks)$density )$p.value Warning message: In chisq.test(hist(pnorm(x), plot = FALSE, breaks = breaks)$density) : Chi-squared approximation may be incorrect > if (pval < alpha) stop("chisq test FAILED! p-value=",signif(pval)) > rm(unr) > > ## use logPDF (use ARS to test 'print' function) > gausspdf <- function (x) { -0.5*x^2 } > gaussdpdf <- function (x) { -x } > gauss <- new("unuran.cont", pdf=gausspdf, dpdf=gaussdpdf, islog=TRUE, lb=-Inf, ub=Inf) > unr <- unuran.new(gauss, "ars") > unr Object is UNU.RAN object: method: ars distr: [S4 class] inversion: FALSE generator ID: ARS.007 distribution: name = unknown type = continuous univariate distribution functions = logPDF dlogPDF domain = (-inf, inf) method: ARS (Adaptive Rejection Sampling) T_c(x) = log(x) ... c = 0 performance characteristics: area(hat) = 4.09235 [ log = 1.40912 ] rejection constant = 1.634 [approx.] # intervals = 3 > x <- unuran.sample(unr, samplesize) > pval <- chisq.test( hist(pnorm(x),plot=FALSE,breaks=breaks)$density )$p.value Warning message: In chisq.test(hist(pnorm(x), plot = FALSE, breaks = breaks)$density) : Chi-squared approximation may be incorrect > if (pval < alpha) stop("chisq test FAILED! p-value=",signif(pval)) > rm(unr) > > > ## --- Discrete distributions ----------------------------------------------- > > ## Create an object > unr <- new("unuran", "binomial(20,0.5)", "dgt") > > ## Draw samples > unuran.sample(unr) [1] 10 > unuran.sample(unr,10) [1] 8 9 8 11 12 8 8 8 7 15 > x <- unuran.sample(unr, samplesize) > rm(unr) > > > ## --- Discrete distributions - S4 distribution object ---------------------- > > ## use PV > pv <- dbinom(0:100,100,0.3) > binom <- new("unuran.discr",pv=pv,lb=0) > unr <- unuran.new(binom, "dgt") > x <- unuran.sample(unr, samplesize) > pval <- chisq.test( hist(pbinom(x,100,0.3),plot=FALSE)$density )$p.value Warning message: In chisq.test(hist(pbinom(x, 100, 0.3), plot = FALSE)$density) : Chi-squared approximation may be incorrect > if (pval < alpha) stop("chisq test FAILED! p-value=",signif(pval)) > rm(unr) > > ## use PMF > pmf <- function(x) dbinom(x,100,0.3) > binom <- new("unuran.discr",pmf=pmf,lb=0,ub=100) > unr <- unuran.new(binom, "dgt") > x <- unuran.sample(unr, samplesize) > pval <- chisq.test( hist(pbinom(x,100,0.3),plot=FALSE)$density )$p.value Warning message: In chisq.test(hist(pbinom(x, 100, 0.3), plot = FALSE)$density) : Chi-squared approximation may be incorrect > if (pval < alpha) stop("chisq test FAILED! p-value=",signif(pval)) > rm(unr) > > ## use PMF > pmf <- function(x) dbinom(x,100,0.3) > binom <- new("unuran.discr",pmf=pmf,lb=0,ub=100) > unr <- unuran.new(binom, "dari") > x <- unuran.sample(unr, samplesize) > pval <- chisq.test( hist(pbinom(x,100,0.3),plot=FALSE)$density )$p.value Warning message: In chisq.test(hist(pbinom(x, 100, 0.3), plot = FALSE)$density) : Chi-squared approximation may be incorrect > if (pval < alpha) stop("chisq test FAILED! p-value=",signif(pval)) > rm(unr) > > > ## --- Continuous Multivariate distributions -------------------------------- > > mvpdf <- function (x) { exp(-sum(x^2)) } > mvd <- new("unuran.cmv", dim=2, pdf=mvpdf, mode=c(0,0)) > unr <- unuran.new(mvd, "hitro") > x <- unuran.sample(unr, 10) > x [,1] [,2] [1,] -0.1413374 0.0000000 [2,] -0.1413374 -0.2489752 [3,] -0.1413505 -0.2489984 [4,] 0.3152893 -0.2489984 [5,] 0.3152893 -0.6237342 [6,] 0.4744190 -0.9385394 [7,] -0.2498318 -0.9385394 [8,] -0.2498318 0.1922014 [9,] -0.2548681 0.1960760 [10,] -0.6784643 0.1960760 > rm(unr) > > unr <- unuran.new(mvd, "vnrou") > x <- unuran.sample(unr, 10) > x [,1] [,2] [1,] 0.4175291 0.0148600 [2,] 1.0662819 -0.4656075 [3,] 0.0814565 -0.4663378 [4,] -0.9823273 1.0288373 [5,] 0.1229493 0.8051288 [6,] 0.4825418 -1.0448065 [7,] 0.4657918 1.6286097 [8,] -1.8204878 -0.1462804 [9,] -0.5039755 0.2405246 [10,] 0.9659430 1.6582299 > rm(unr) > > > ## --- quantile function ---------------------------------------------------- > > ## test U-error > unr <- unuran.new("normal()","hinv; u_resolution=1.e-13") > > ## single double as argument > Tmax <- 0 > for (U in (0:20)/20) { + T <- pnorm( uq(unr,U) ) - U + print (T) + Tmax <- max(abs(T),Tmax) + } [1] 0 [1] 3.985701e-14 [1] 2.378653e-14 [1] 1.898481e-14 [1] 2.581269e-15 [1] 1.49325e-14 [1] 1.29341e-14 [1] 9.214851e-15 [1] 6.383782e-15 [1] 4.518608e-14 [1] -5.551115e-17 [1] -4.52971e-14 [1] -6.439294e-15 [1] -9.214851e-15 [1] -1.298961e-14 [1] -1.509903e-14 [1] -2.664535e-15 [1] -1.898481e-14 [1] -2.38698e-14 [1] -3.985701e-14 [1] 0 > cat("Max. error =",Tmax,"\n") Max. error = 4.52971e-14 > if (Tmax > 1.e-13) stop ("Max. error exceeds limit") > > ## vector argument > U <- (0:20)/20 > T <- pnorm( uq(unr,U) ) - U > print (T) [1] 0.000000e+00 3.985701e-14 2.378653e-14 1.898481e-14 2.581269e-15 [6] 1.493250e-14 1.293410e-14 9.214851e-15 6.383782e-15 4.518608e-14 [11] -5.551115e-17 -4.529710e-14 -6.439294e-15 -9.214851e-15 -1.298961e-14 [16] -1.509903e-14 -2.664535e-15 -1.898481e-14 -2.386980e-14 -3.985701e-14 [21] 0.000000e+00 > Tmax <- max(abs(T)) > cat("Max. error =",Tmax,"\n") Max. error = 4.52971e-14 > if (Tmax > 1.e-13) stop ("Max. error exceeds limit") > > ## special arguments > for (U in c(-1,-0.001,0,0.5,1.,1.001,NA,NaN) ) { + cat ("U =",U,"\tX =",uq(unr,U),"\n") + } [UNU.RAN - warning] argument out of domain: U not in [0,1] U = -1 X = -Inf [UNU.RAN - warning] argument out of domain: U not in [0,1] U = -0.001 X = -Inf U = 0 X = -Inf U = 0.5 X = -1.387779e-16 U = 1 X = Inf [UNU.RAN - warning] argument out of domain: U not in [0,1] U = 1.001 X = Inf U = NA X = NA U = NaN X = NaN > > U <- c(-1,-0.001,0,0.5,1.,1.001,NA,NaN) > T <- uq(unr,U) [UNU.RAN - warning] argument out of domain: U not in [0,1] [UNU.RAN - warning] argument out of domain: U not in [0,1] [UNU.RAN - warning] argument out of domain: U not in [0,1] > rbind(U,T) [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] U -1 -0.001 0 5.000000e-01 1 1.001 NA NaN T -Inf -Inf -Inf -1.387779e-16 Inf Inf NA NaN > > uq(unr,numeric()) numeric(0) > > rm(unr) > > ## test whether 'uq' throws an error when UNU.RAN object does not implement > ## an inversion method > unr <- unuran.new("normal()","tdr") > if( ! is.error( uq(unr,0.5) ) ) + stop("'uq' does not detect UNU.RAN object with non-inversion method") Error in uq(unr, 0.5) : [UNU.RAN - error] invalid UNU.RAN object: inversion method required! Use methods 'HINV', 'NINV', 'PINV'; or 'DGT' > > ## test whether 'uq' detects invalid arguments > if( ! is.error( uq(1,0.5) ) ) + stop("'uq' does not detect invalid argument 'unr'") Error in uq(1, 0.5) : argument 'unr' must be UNU.RAN object; method not implemented for distribution objects > > if( ! is.error( uq(unr,"a") ) ) + stop("'uq' does not detect invalid argument 'U'") Error in uq(unr, "a") : [UNU.RAN - error] argument invalid: 'U' must be numeric > > rm(unr) > > > ## --- density function ----------------------------------------------------- > > distr <- unuran.cont.new(pdf=function(x){exp(-x)}, lb=0,ub=Inf) > x <- rexp(100) > e <- max(abs(ud(distr, x) - dexp(x))) > e; if (e>1.e-10) stop("error too large") [1] 0 > > distr <- unuran.cont.new(pdf=function(x){-x}, islog=TRUE, lb=0,ub=Inf) > x <- rexp(100) > e <- max(abs(ud(distr, x, islog=TRUE) + x)) > e; if (e>1.e-10) stop("error too large") [1] 0 > > distr <- udgeom(prob=0.8) > x <- rgeom(100, prob=0.8) > e <- max(abs(ud(distr,x) - dgeom(x,prob=0.8))) > e; if (e>1.e-10) stop("error too large") [1] 0 > > rm(distr,x,e) > > unr <- pinv.new(pdf=function(x){exp(-x)}, lb=0,ub=Inf) > x <- rexp(100) > e <- max(abs(ud(unr, x) - dexp(x))) > e; if (e>1.e-10) stop("error too large") [1] 0 > > unr <- pinv.new(pdf=function(x){-x}, islog=TRUE, lb=0,ub=Inf) > x <- rexp(100) > e <- max(abs(ud(unr, x, islog=TRUE) +x)) > e; if (e>1.e-10) stop("error too large") [1] 0 > > unr <- darid.new(udgeom(prob=0.8)) > x <- rgeom(100, prob=0.8) > e <- max(abs(ud(unr,x) - dgeom(x,prob=0.8))) > e; if (e>1.e-10) stop("error too large") [1] 0 > > rm(unr,x,e) > > > distr <- unuran.cont.new(lb=0,ub=1) > if (!all(is.na(ud(distr,1)))) + stop("'ud' ignores missing PDF") Warning message: In ud(distr, 1) : [UNU.RAN - error] UNU.RAN object does not contain (log)PDF > > distr <- unuran.discr.new(lb=0,ub=1) > if (!all(is.na(ud(distr,1)))) + stop("'ud' ignores missing PMF") Warning message: In ud(distr, 1) : [UNU.RAN - error] UNU.RAN object does not contain (log)PMF > > unr <- pinv.new(cdf=pexp,lb=0,ub=Inf) > if (!all(is.na(ud(distr,1)))) + stop("'ud' ignores missing PDF") Warning message: In ud(distr, 1) : [UNU.RAN - error] UNU.RAN object does not contain (log)PMF > > unr <- pinv.new(pdf=dexp,lb=0,ub=Inf) > unuran.packed(unr) <- TRUE > if( ! is.error( ud(unr, 1) ) ) + stop("'ud' does not detect packed generator object") Error in ud(unr, 1) : [UNU.RAN - error] cannot compute PDF for packed UNU.RAN object > > rm(distr,unr) > > > ## --- distribution function ------------------------------------------------ > > distr <- unuran.cont.new(cdf=function(x){1-exp(-x)}, lb=0,ub=Inf) > x <- rexp(100) > e <- max(abs(up(distr, x) - pexp(x))) > e; if (e>1.e-10) stop("error too large") [1] 5.551115e-17 > > distr <- udgeom(prob=0.8) > x <- rgeom(100, prob=0.8) > e <- max(abs(up(distr,x) - pgeom(x,prob=0.8))) > e; if (e>1.e-10) stop("error too large") [1] 0 > > unr <- pinv.new(pdf=function(x){exp(5-x)}, lb=0,ub=Inf) > x <- rexp(100) > e <- max(abs(up(unr, x) - pexp(x))) > e; if (e>1.e-10) stop("error too large") [1] 4.368395e-12 > > rm(distr,unr,x,e) > > > distr <- unuran.cont.new(lb=0,ub=1) > if( ! is.error( up(distr,1) ) ) + stop("'up' ignores missing CDF") Error in up(distr, 1) : [UNU.RAN - error] UNU.RAN object does not contain CDF > > distr <- unuran.discr.new(lb=0,ub=1) > if( ! is.error( up(distr,1) ) ) + stop("'up' ignores missing CDF") Error in up(distr, 1) : [UNU.RAN - error] UNU.RAN object does not contain CDF > > unr <- tdr.new(pdf=dexp,lb=0,ub=Inf) > if( ! is.error( up(unr,1) ) ) + stop("'up' ignores invalid method PINV") Error in up(unr, 1) : [UNU.RAN - error] function requires method PINV > > unr <- pinv.new(pdf=dexp,lb=0,ub=Inf) > unuran.packed(unr) <- TRUE > if( ! is.error( up(unr, 1) ) ) + stop("'up' does not detect packed generator object") Error in up(unr, 1) : [UNU.RAN - error] cannot compute CDF for packed UNU.RAN object > > rm(distr,unr) > > > ## --- pack ----------------------------------------------------------------- > > ## check print with unuran.packed > unr <- pinv.new(dnorm,lb=0,ub=Inf) > unuran.packed(unr) <- TRUE > unr Object is UNU.RAN object: method: pinv;usepdf;u_resolution=1e-10;smoothness=0;keepcdf=on distr: [S4 class] inversion: TRUE Object is PACKED ! > unuran.details(unr) Object is UNU.RAN object: method: pinv;usepdf;u_resolution=1e-10;smoothness=0;keepcdf=on distr: [S4 class] inversion: TRUE Object is PACKED ! > unuran.details(unr,show=TRUE,return.list=TRUE) Object is UNU.RAN object: method: pinv;usepdf;u_resolution=1e-10;smoothness=0;keepcdf=on distr: [S4 class] inversion: TRUE Object is PACKED ! Object is PACKED ! > unuran.details(unr,show=FALSE,return.list=TRUE) Object is PACKED ! > print(unuran.details(unr,show=FALSE,return.list=TRUE)) Object is PACKED ! NULL > print(unuran.details(unr,show=FALSE,debug=TRUE)) Object is PACKED ! NULL > unuran.details(unr,show=FALSE,return.list=FALSE) > rm(unr) > > ## check whether un/packing works/fails > unr <- pinv.new(dnorm,lb=0,ub=Inf) > ## this should be o.k. > unuran.packed(unr) <- FALSE > ## this should fail with error > unuran.packed(unr) <- TRUE > if( ! is.error( unuran.packed(unr) <- FALSE ) ) + stop("'unuran.packed' tries to unpack UNU.RAN objects") Error : [UNU.RAN - error] Cannot unpack 'unuran' object > ## this should be o.k. > unuran.packed(unr) <- TRUE Warning message: [UNU.RAN - warning] object already PACKED > rm(unr) > > ## test whether non-packable objects are treated correctly > unr <- tdr.new(dnorm,lb=0,ub=Inf) > if( ! is.error( unuran.packed(unr) <- TRUE ) ) + stop("'unuran.packed' does not detect non-packable UNU.RAN objects") Error : [UNU.RAN - error] cannot pack UNU.RAN object > rm(unr) > > > ## --- mixture -------------------------------------------------------------- > > comp <- c(unuran.new("normal"),unuran.new("cauchy"),unuran.new("exponential")) > prob <- c(1,2,3) > unr <- mixt.new(prob,comp) > x <- unuran.sample(unr, 10) > x [1] 1.7645417 2.1617851 0.8693685 1.4906015 -0.1282619 1.3519479 [7] 2.6755241 -0.6539987 -2.6498065 0.8089170 > rm(unr) > > comp <- c(pinvd.new(udnorm(lb=-Inf,ub=-1)), + pinvd.new(udcauchy(lb=-1,ub=1)), + pinvd.new(udexp(lb=1,ub=Inf)) ) > prob <- c(1,2,3) > unr <- mixt.new(prob,comp,inversion=TRUE) > x <- unuran.sample(unr, 10) > x [1] -1.07237310 1.82810164 0.01126635 -1.45231972 1.76962315 1.04170946 [7] 1.34244672 -0.11644349 1.17637391 2.73140543 > rm(unr) > > > ## --- Check generator ------------------------------------------------------ > > ## run with defaults > unr <- tdrd.new(udnorm()) > unuran.verify.hat(unr) Check inequality squeeze(x) <= density(x) <= hat(x) for automatic rejection method: 0 out of 1e+05 (= 0%) points failed. No problems have been detected! > rm(unr) > > ## do not show result > unr <- tdrd.new(udnorm()) > unuran.verify.hat(unr,show=FALSE) > rm(unr) > > ## another example > pdf <- function(x) { -x - 0.999*log(x) } > dpdf <- function(x) { -1 - 0.999/x } > distr <- unuran.cont.new(pdf=pdf, dpdf=dpdf, islog=TRUE, lb=0, ub=0.5, mode=0) > unr <- unuran.new(distr,method="itdr;cp=-0.999") > unuran.verify.hat(unr) Check inequality squeeze(x) <= density(x) <= hat(x) for automatic rejection method: 0 out of 1e+05 (= 0%) points failed. No problems have been detected! > rm(unr); rm(distr); rm(pdf); rm(dpdf) > > ## example where the hat violates condition hat(x) >= pdf(x) > ## it also stores the ratio in variable 'failed' > unr <- unuran.new(udnorm(),method="nrou;v=0.6") > failed <- unuran.verify.hat(unr) Check inequality squeeze(x) <= density(x) <= hat(x) for automatic rejection method: 33445 out of 1e+05 (= 33.44%) points failed. Verifying hat and squeeze failed! The generator does not sample from the requested distribution! > failed [1] 0.33445 > rm(unr); rm (failed) > > ## Example for a method that does not implement rejection > unr <- pinvd.new(udnorm()) > if (! is.error( unuran.verify.hat(unr) ) ) + stop("method must throw error for inversion method") Error : [UNU.RAN - error] Method not supported > rm(unr) > > > ## --- End ------------------------------------------------------------------ > > detach("package:Runuran",unload = TRUE) > > ## -------------------------------------------------------------------------- > > proc.time() user system elapsed 6.10 0.43 6.53