[
https://issues.apache.org/jira/browse/RNG-195?page=com.atlassian.jira.plugin.system.issuetabpanels:comment-tabpanel&focusedCommentId=18126326#comment-18126326
]
Alex Herbert commented on RNG-195:
----------------------------------
Examples of output from the updated Dirichlet sampler. The mean and variance
have been computed using the moments of the distribution ([Dirichlet moments
(Wikipedia)|https://en.wikipedia.org/wiki/Dirichlet_distribution#Moments]).
This was done using BigDecimal to support extreme alpha values without
over/underflow.
h2. Symmetric small alpha
When alpha is the same across categories, the sample will concentrate towards 1
category as alpha reduces. When alpha is small (<1e-5) the samples appear to be
a random choice of one of the categories. I have verified that across many
samples the frequency of observed categories is uniform using a chi-square
test. At extreme concentrations the variance is similar to small concentrations
and variation is still observed in the sample.
h3. alpha: [0.1, 0.1, 0.1]
{noformat}
m: 0.3333333333333333 0.3333333333333333 0.3333333333333333
v: 0.1709401709401709 0.1709401709401709 0.1709401709401709
[0.15809228819213816, 5.636711602005743E-6, 0.8419020750962597]
[3.394144793503779E-8, 0.043771340315133085, 0.956228625743419]
[2.1360623407495405E-4, 0.9447960133996312, 0.05499038036629393]
[1.7158145387496422E-4, 0.9527804893891213, 0.047047929157003864]
[1.391774435867783E-4, 2.7499931349608023E-5, 0.9998333226250636]
[3.9772778487113543E-10, 1.7264176742667857E-5, 0.9999827354255295]
[0.021171906599323942, 0.045555499706225785, 0.9332725936944503]
[7.652207114103051E-6, 0.9982696361270899, 0.001722711665796009]
[0.32495911046589476, 0.014534884122104514, 0.6605060054120008]
[1.8602097612115832E-11, 0.5486673208321027, 0.4513326791492952]
{noformat}
h3. alpha: [0.01, 0.01, 0.01]
{noformat}
m: 0.3333333333333333 0.3333333333333333 0.3333333333333333
v: 0.2157497303128371 0.2157497303128371 0.2157497303128371
[0.0020737720573923042, 1.4092012907875255E-32, 0.9979262279426077]
[2.1658726101085867E-154, 7.453619917873339E-7, 0.9999992546380081]
[4.428656424244222E-4, 0.9995316277924401, 2.5506565135465564E-5]
[3.514815283932834E-65, 0.9999999993124695, 6.875305174293745E-10]
[4.881412756375219E-40, 1.062003410544296E-75, 1.0]
[0.9999999999999274, 2.0374749383977687E-49, 7.264906031144927E-14]
[0.9592889825470766, 0.04071101745292331, 5.659907229325881E-48]
[0.9682831394890384, 1.4201725139846832E-53, 0.03171686051096167]
[3.3783638198645964E-32, 0.43402445982976123, 0.5659755401702388]
[3.6148356233351014E-8, 7.916494795361955E-30, 0.9999999638516437]
{noformat}
h3. alpha: [0.001, 0.001, 0.001]
{noformat}
m: 0.3333333333333333 0.3333333333333333 0.3333333333333333
v: 0.2215575495735017 0.2215575495735017 0.2215575495735017
[7.0559501300893855E-276, 1.0, 2.4955696190414105E-84]
[1.0, 0.0, 0.0]
[0.0, 1.0, 9.484491006193E-243]
[1.6141250965640683E-81, 3.0926370203279035E-260, 1.0]
[0.0, 0.0, 1.0]
[9.521045197744373E-195, 3.981843972245724E-235, 1.0]
[1.0, 0.0, 4.202778794986451E-204]
[0.9999999999999993, 5.939030291354542E-16, 1.1504289134429875E-66]
[0.0, 0.0, 1.0]
[1.0, 0.0, 9.52355287715784E-33]
{noformat}
h3. alpha: [1.0E-5, 1.0E-5, 1.0E-5]
{noformat}
m: 0.3333333333333333 0.3333333333333333 0.3333333333333333
v: 0.2222155557555496 0.2222155557555496 0.2222155557555496
[0.0, 0.0, 1.0]
[1.0, 0.0, 0.0]
[0.0, 0.0, 1.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[0.0, 0.0, 1.0]
[0.0, 0.0, 1.0]
[0.0, 1.0, 0.0]
[0.0, 1.0, 0.0]
[1.0, 0.0, 0.0]
{noformat}
h3. alpha: [4.9E-324, 4.9E-324, 4.9E-324]
{noformat}
m: 0.3333333333333333 0.3333333333333333 0.3333333333333333
v: 0.2222222222222222 0.2222222222222222 0.2222222222222222
[0.0, 1.0, 0.0]
[0.0, 0.0, 1.0]
[0.0, 1.0, 0.0]
[1.0, 0.0, 0.0]
[0.0, 1.0, 0.0]
[0.0, 0.0, 1.0]
[0.0, 0.0, 1.0]
[0.0, 0.0, 1.0]
[0.0, 0.0, 1.0]
[1.0, 0.0, 0.0]
{noformat}
h2. Asymmetric small alpha
The distribution is biased towards the category with the largest alpha. With a
10-fold difference of concentrations one category dominates the samples. At
extreme concentrations the variance is similar to small concentrations and
variation is still observed in the sample.
h3. alpha: [0.1, 0.01, 0.001]
{noformat}
m: 0.9009009009009009 0.09009009009009009 0.009009009009009009
v: 0.08035865675683711 0.07378385756764134 0.008035865675683712
[0.8946038279389222, 0.10539617206107796, 4.200758506266063E-93]
[1.0, 4.298107611169935E-25, 0.0]
[1.0, 5.067756010693927E-57, 0.0]
[1.0, 2.2889648826401724E-71, 0.0]
[1.0, 3.8604434266769892E-56, 4.1752765817728853E-196]
[0.9999999999999999, 1.0959018489777459E-16, 1.47170303879052E-78]
[1.0, 8.763719536495225E-39, 0.0]
[0.9999999996975344, 1.1228209890802343E-17, 3.0246550569607727E-10]
[0.9999999999997584, 2.4171618690469716E-13, 0.0]
[0.9999999948795391, 5.120460979564197E-9, 3.5105824922422084E-210]
{noformat}
h3. alpha: [1.0E-5, 1.0E-6, 1.0E-7]
{noformat}
m: 0.9009009009009009 0.09009009009009009 0.009009009009009007
v: 0.08927747667685491 0.08197295585783951 0.008927747667685492
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[0.0, 1.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[0.0, 1.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
{noformat}
h3. alpha: [1.0E-301, 1.0E-302, 1.0E-303]
{noformat}
m: 0.9009009009009009 0.09009009009009009 0.009009009009009007
v: 0.08927846765684604 0.08197386575764953 0.008927846765684602
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[0.0, 1.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
[1.0, 0.0, 0.0]
{noformat}
h2. Symmetric large alpha
Samples increasingly tend towards 1/k for all k categories when alpha
increases. At large alpha the variance is sub-normal and the sample is uniform.
Note the sampler outputs uniform samples from around alpha 1e33 due to limits
of floating point arithmetic.
h3. alpha: [1.0E10, 1.0E10, 1.0E10]
{noformat}
m: 0.3333333333333333 0.3333333333333333 0.3333333333333333
v: 7.407407407160493E-12 7.407407407160493E-12 7.407407407160493E-12
[0.3333342434133697, 0.3333307200299129, 0.3333350365567174]
[0.3333304564571932, 0.333336813988147, 0.33333272955465976]
[0.333335725193178, 0.3333349765641703, 0.33332929824265173]
[0.3333373707285592, 0.33333368515905404, 0.33332894411238684]
[0.3333322364354621, 0.3333363771476584, 0.3333313864168794]
[0.33333063103778504, 0.33333756859514097, 0.33333180036707405]
[0.3333378288103135, 0.3333282375475025, 0.33333393364218394]
[0.33333001288095454, 0.333337003947986, 0.33333298317105947]
[0.3333331108399213, 0.33333817170277313, 0.33332871745730563]
[0.3333286726477712, 0.3333363734702074, 0.3333349538820214]
{noformat}
h3. alpha: [1.0E20, 1.0E20, 1.0E20]
{noformat}
m: 0.3333333333333333 0.3333333333333333 0.3333333333333333
v: 7.407407407407407E-22 7.407407407407407E-22 7.407407407407407E-22
[0.33333333334410514, 0.33333333329843623, 0.3333333333574586]
[0.33333333334117454, 0.3333333333185414, 0.33333333334028403]
[0.33333333328826126, 0.33333333333837106, 0.3333333333733677]
[0.33333333327935655, 0.33333333339785376, 0.3333333333227897]
[0.3333333333761238, 0.33333333328492326, 0.33333333333895293]
[0.3333333333222844, 0.33333333336602067, 0.33333333331169496]
[0.3333333333210733, 0.3333333333734735, 0.3333333333054532]
[0.3333333333618135, 0.33333333334510623, 0.3333333332930803]
[0.3333333332974264, 0.3333333333454568, 0.3333333333571168]
[0.33333333330313686, 0.3333333333546322, 0.3333333333422309]
{noformat}
h3. alpha: [1.0E25, 1.0E25, 1.0E25]
{noformat}
m: 0.3333333333333333 0.3333333333333333 0.3333333333333333
v: 7.407407407407407E-27 7.407407407407407E-27 7.407407407407407E-27
[0.33333333333331283, 0.3333333333333697, 0.33333333333331755]
[0.3333333333333586, 0.3333333333331833, 0.33333333333345805]
[0.3333333333334557, 0.3333333333332757, 0.3333333333332686]
[0.3333333333333278, 0.3333333333333278, 0.33333333333334436]
[0.33333333333333964, 0.33333333333322834, 0.333333333333432]
[0.3333333333333894, 0.3333333333334415, 0.3333333333331691]
[0.33333333333334203, 0.3333333333332899, 0.33333333333336806]
[0.33333333333334836, 0.33333333333343595, 0.3333333333332157]
[0.3333333333332583, 0.333333333333282, 0.33333333333345966]
[0.33333333333324094, 0.3333333333334352, 0.3333333333333239]
{noformat}
h3. alpha: [1.7976931348623157E308, 1.7976931348623157E308,
1.7976931348623157E308]
{noformat}
m: 0.3333333333333333 0.3333333333333333 0.3333333333333333
v: 4.1205071453837E-310 4.1205071453837E-310 4.1205071453837E-310
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
[0.3333333333333333, 0.3333333333333333, 0.3333333333333333]
{noformat}
h2. Asymmetric large alpha
The largest alpha dominates the sample. When all alphas are above 1e33 the
sample is fixed. The sample value in each category is alpha normalized by the
sum of alphas. This can be computed using the log-sum-exp method. It varies
from the expected mean by a few units in the last place (ULP).
h3. alpha: [100000.0, 1000000.0, 1.0E7]
{noformat}
m: 0.009009009009009009 0.09009009009009009 0.9009009009009009
v: 8.043104469706447E-10 7.385032285821374E-9 8.043104469706447E-9
[0.009025984974442375, 0.08979819497593673, 0.9011758200496208]
[0.009028777740323358, 0.09007653118656186, 0.9008946910731148]
[0.009012366023404564, 0.0902312806035209, 0.9007563533730745]
[0.009004002623982819, 0.09037626951813134, 0.9006197278578857]
[0.008989470846643513, 0.0901601112436722, 0.9008504179096842]
[0.008976376781218312, 0.09009522094872002, 0.9009284022700615]
[0.008960715676146903, 0.0899869029246212, 0.9010523813992318]
[0.008986046333862194, 0.09013932499544734, 0.9008746286706906]
[0.008974126496832137, 0.09007996306170923, 0.9009459104414587]
[0.009015324051070433, 0.09005582764640416, 0.9009288483025254]
{noformat}
h3. alpha: [1.0E15, 1.0E16, 1.0E17]
{noformat}
m: 0.009009009009009009 0.09009009009009009 0.9009009009009009
v: 8.043105194310454E-20 7.385032951139598E-19 8.043105194310454E-19
[0.0090090089995922, 0.09009008908481204, 0.9009009019155958]
[0.00900900886528343, 0.09009009093870446, 0.900900900196012]
[0.009009009292973028, 0.0900900903401551, 0.9009009003668719]
[0.009009008637200926, 0.09009008882540358, 0.9009009025373955]
[0.009009009741716135, 0.09009009048977117, 0.9009008997685126]
[0.009009009217773032, 0.09009008886948855, 0.9009009019127385]
[0.009009009097927081, 0.09009008878154719, 0.9009009021205256]
[0.00900900935498803, 0.09009008927036177, 0.9009009013746502]
[0.009009009148904103, 0.09009009003757866, 0.9009009008135171]
[0.009009008851155603, 0.09009008898974978, 0.9009009021590947]
{noformat}
h3. alpha: [1.0E25, 1.0E26, 1.0E27]
{noformat}
m: 0.00900900900900901 0.09009009009009009 0.9009009009009009
v: 8.043105194310454E-30 7.385032951139598E-29 8.043105194310454E-29
[0.00900900900900325, 0.09009009009008732, 0.9009009009009095]
[0.009009009009008589, 0.09009009009009847, 0.9009009009008929]
[0.00900900900901231, 0.09009009009010048, 0.9009009009008873]
[0.009009009009012437, 0.09009009009008126, 0.9009009009009064]
[0.009009009009005798, 0.09009009009007952, 0.9009009009009147]
[0.009009009009009847, 0.09009009009009313, 0.900900900900897]
[0.009009009009009333, 0.09009009009008735, 0.9009009009009034]
[0.009009009009011288, 0.09009009009010754, 0.9009009009008812]
[0.009009009009009984, 0.09009009009009834, 0.9009009009008916]
[0.009009009009006601, 0.09009009009008177, 0.9009009009009117]
{noformat}
h3. alpha: [1.0E35, 1.0E36, 1.0E37]
{noformat}
m: 0.009009009009009009 0.0900900900900901 0.9009009009009009
v: 8.043105194310454E-40 7.3850329511396E-39 8.043105194310454E-39
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
[0.009009009009009056, 0.0900900900900897, 0.9009009009009011]
{noformat}
> DirichletSampler has unbounded recursion and stack overflow for small alpha
> parameters
> --------------------------------------------------------------------------------------
>
> Key: RNG-195
> URL: https://issues.apache.org/jira/browse/RNG-195
> Project: Commons RNG
> Issue Type: Bug
> Components: sampling
> Affects Versions: 1.3
> Reporter: Alex Herbert
> Priority: Minor
>
> The DirichletSampler uses a rejection method to sample. A gamma sample is
> obtained for each concentration parameter alpha. The vector of k samples is
> normalised to unit length. When the sum of the samples is zero or infinite
> the sample cannot be normalised and recursion occurs to generate another
> sample.
> If the alpha parameters are very small the gamma samples are always zero and
> a stack overflow error occurs.
> Note that recursion is expected to be used to indicate a non-functional RNG
> underlying the sampler. It was not the original intention to use a stack
> overflow error to indicate a badly parameterised sampler where samples are
> impossible.
> Issue identified using a security scan.
--
This message was sent by Atlassian Jira
(v8.20.10#820010)
