Measuring Kernel Density

As previously discussed, we can compute the expected kernel density S0 of a mask list of size Q which is usually but not necessarily a power of 2. But this is only useful if we can measure the actual kernel density S of an arbitrary such list.

To do this, let's revisit RandomMapGet[] and declare it at the "In[1]:=" prompt:
__________________________________________________________________________
In[1]:= RandomMapGet[bitcount_]:=(intcount0=BitShiftLeft[1, bitcount]; intmax0=intcount0-1; randommap0=Array[RandomInteger[intmax0]&, intcount0]; Return[randommap0])

In[2]:= bitcount=16;

In[3]:= randommap=RandomMapGet[bitcount];
__________________________________________________________________________

Now randommap contains (Q = (2 ^ 16)) 16-bit masks.

How large it its kernel? To answer this question, we need to discover which of all possible 16-bit masks come back to themselves, iteratively speaking. In other words, how many of them are on limit cycles?

We could in theory try every single mask and test whether or not it returns to itself in at most Q iterations. But that would take a long time.

A faster way to do this is to use double buffers, notionally map X and map Y. Begin by creating bitmap Z consisting of all zeroes. Next set map X to the input map, which is randommap in this case. Start of loop: Now, one index at a time, set Y[i] to X[X[i]] (assuming zero-based indexes, which is annoyingly not the case in Mathematica). After Y has been completely populated thus, do the reverse: one index at a time, set X[i] to Y[Y[i]]. Now copy Z to Z'. Then evaluate the "profile" of X: set bits in Z corresponding to masks present in X; if the mask is present one or more times, then the bit is one, else zero. Loop back to "start of loop" until (Z = Z'). Count the number of ones in Z and divide by Q; this is the kernel density S.

This could probably be optimized further, but it's fast enough for some useful minimodel analysis. Here's the code:
__________________________________________________________________________

KernelDensityGet[bitcount_, samplemap_]:=(intcount1=BitShiftLeft[1, bitcount]; iteratedmapa1=samplemap; iteratedmapb1=samplemap; bitmap1=BitShiftLeft[1, intcount1]; profilebitmap1=bitmap1; profilebitmapprev1=BitOr[bitmap1, 1]; While[profilebitmap1!=profilebitmapprev1, profilebitmapprev1=profilebitmap1; For[i1=1, i1<=intcount1, i1++, iteratedmapb1[[i1]]=iteratedmapa1[[iteratedmapa1[[i1]]+1]]]; For[i1=1, i1<=intcount1, i1++, iteratedmapa1[[i1]]=iteratedmapb1[[iteratedmapb1[[i1]]+1]]]; profilebitmap1=bitmap1; For[i1=1, i1<=intcount1, i1++, bit1=iteratedmapa1[[i1]]; If[BitGet[profilebitmap1, bit1]==0, profilebitmap1=BitSet[profilebitmap1, bit1]]]]; s1=(DigitCount[profilebitmap1, 2, 1]-1)/intcount1+0.; Return[s1])

__________________________________________________________________________

Now run it on the previously generated randommap (your answer will vary):

__________________________________________________________________________

In[4]:= s=KernelDensityGet[bitcount, randommap]

Out[4]= 0.0055542
__________________________________________________________________________

You can compute S0 from KernelDensityExpectedGet[] at the aforelinked post. It turns out to be 0.00489068 for 16 bits. Thus the actual kernel density in my one random test above exceeded expectation by about 14%. Unlike the dyspoissonism floor, which is a lower bound, actual kernel density could be greater or lesser than the expected value. I have no idea whether or not 14% is "good enough", but the important thing is that (1) by iterating this process many times, we could obtain a very accurate estimate of adherence to expected kernel density and (2) we have yet another "hard" randomness evaluation tool in our hands by which to measure the relative merits of various PRNGs. Furthermore, kernel skew provides a normalized method of comparing the kernel density among disparate mask lists.