(*
 * Implementation of NIST Special Publication 800-22
 * "Statistical Test Suite for Random and Psuedornadom
 *  Number Generators for Cryptographic Applications"
 * 
 * http://modp.com/release/nist800-22/
 * Version 0.9.2 - 07-Oct-2006
 * Nick Galbreath  nickg [at] modp [dot] com
 * Copyright 2005, 2006
 * All rights reserved.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions are
 * met:
 *
 *   Redistributions of source code must retain the above copyright
 *   notice, this list of conditions and the following disclaimer.
 *
 *   Redistributions in binary form must reproduce the above copyright
 *   notice, this list of conditions and the following disclaimer in the
 *   documentation and/or other materials provided with the distribution.
 *
 *   Neither the name of the modp.com nor the names of its
 *   contributors may be used to endorse or promote products derived from
 *   this software without specific prior written permission.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 * OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 *
 * This is the standard "new" BSD license:
 * http://www.opensource.org/licenses/bsd-license.php
 *)

WARNING

This is a preview release.  This needs to be turned into a package and documented. Interfaces and names will change. 07-Oct-2006

Common Code

<< Statistics`DataManipulation`
\(\(\(<< Statistics`NormalDistribution`\)\(
\)\(
\)
    \) (*\ Using\ RegExp\ is\ fastest, \ 
      but\ can\ only\ handle\ limited\ size\ *) \)
\(\(\(StringPartition[s_String, \ m_Integer]\  := \ 
If[
        m\  < \ 100, 
StringCases[s, \ 
          Apply[StringExpression, \ Table[_, \ {m}]], \ 
          Overlaps -> False], 
FromCharacterCode[
          Partition[ToCharacterCode[s], \ m]]]\)\(
\)
    \)\)
\(StringPartition[s_String, \ m_Integer, \ 
        1]\  := 
StringCases[s, \ 
        Apply[StringExpression, \ Table[_, \ {m}]], \ 
        Overlaps -> True]
\)
BinaryString[x_Integer,  n_Integer]  := 
      StringJoin[
        Map[FromCharacterCode[# + 48]  &, 
          IntegerDigits[x, 2, 
            n]]]

StringPartitionArray[s_String,  m_Integer]  :=  
      FromCharacterCode[
        Partition[ToCharacterCode[s],  m]]

\(RandomBinaryString[m_Integer]\  := \ 
      FromCharacterCode[
        Table[Random[Integer] + 48, \ {m}]]

 (*\ 
      Converts\ a\ string\ of\ "1"\ and\ "0"\ to\ a\ list\ of\ 0\ and\ 1\ *) \)
BinaryStringToDigits[s_String]  :=  
      ToCharacterCode[s]  -  48
BinaryStringToInteger[s_String]  :=  
      FromDigits[ ToCharacterCode[s]  -  48,  2]

2.1 Monobit Frequency Test

\(\(\( (*2.1*) \)\(
\)\(MonobitFrequencyTest[s_String]\  := \ 
      Module[{n\  = \ StringLength[s]}, 
sn\  = 
          Apply[Plus, \((ToCharacterCode[s]\  - \ 48)\)*2 - 
              1]; 
sobs\  = \ 
          Abs[sn]\ /\ \ Sqrt[n]; 
pval = \ 
          Erfc[sobs\ /\ Sqrt[2]]; 
pval\  ≥ \ 
          0.01
]\)\)\)

2.2 Block Frequency Test

\(\(\( (*\ 2.2\ *) \)\(
\)\( (*\ 
      what\ is\ this?\ \ I\ \(forgot!\)\ \ nickg*) \)\(
\)\(BlockFrequencyTest2[
        s_String, \ m_Integer]\  := \ 
      Map[Count[#, 1] &, \ 
        Partition[\((ToCharacterCode[s]\  - \ 48)\), \ 
          m]]


    BlockFrequencyTest[s_String, \ m_Integer]\  := \ 
      Module[{}, 
freq\  = 
          Map[StringCount[#, \ "1"] &, \ StringPartition[s, \ m]]\ /\ 
            m; 
n\  = \ 
          Length[freq]; 
chisqr\  = \ 
          4*m*Apply[Plus, \ 
              Map[\((#\  - \ 1/2)\)^2 &, \ 
                freq]]; 
pval\  = \ 
          GammaRegularized[\ n/2, \ 
            chisqr\ /2]; 
pval\  ≥ \ 
          0.01
]\)\)\)

2.3 Runs Test

\(\(\( (*2.3*) \)\(
\)\(RunsTest[
          s_String]\  := \ 
Module[{n\  = \ 
              StringLength[s]}, 
pi\  = \ 
            StringCount[s, \ "1"]\ /\ n; 
 (*\ 
            the\ regexp\ counts\ number\ of\ runs\ of\ 1  s\ and\ 0  
              s\ *) 
vobs\  = \ 
            StringCount[s, 
              RegularExpression["0+|1+"]]; 
pval\  = \ 
            Erfc[\((vobs\  - \ 2\ *\ n\ *\ pi*\((1 - pi)\))\)/\ \((2*
                    pi*\((1 - pi)\)*
                    Sqrt[2*n])\)]; 
pval\  ≥ \ 0.01
];\)\)\)

2.4 Longest Runs Test

\(\(\( (*2.4*) \)\(
\)\(LongestRunOnes8[
        s_String]\  := \ 
Module[{m = 8, \ k = 3, \ 
          n = \ 16, \ 
          pik\  = \ {0.2148, \ 0.3672, \ 0.2305, \ 
              0.1875}}, 
blocks\  = \ 
          StringPartition[s, \ m]; 
counts\  = \ 
          Map[\ Max[1, \ 
                Min[4, Max[StringLength[StringSplit[#, \ "0"]]]]] &, \ 
            blocks]; 
freqs\  = \ 
          Table[Count[counts, \ i], \ {i, 1, 
              4}]; 
chisqr\  = \ 
          Apply[Plus, 
            MapThread[\ \((\((#1\  - \ 
                          n*\ #2)\)\ ^2)\)\ /\ \((n\ *\ #2)\) &, \ {freqs, \ 
                pik}]]; 
pval\  = \ 
          GammaRegularized[k/2, \ 
            chisqr\ /\ 2]; 
pval\  ≥ \ 
          0.01
]


    LongestRunsOnes128[
        s_String]\  := 
Module[{blocks, m = 128, \ k = 5, \ 
          n = \ 49, \ 
          pik\  = \ {0.1174, \ 0.2430, \ 0.2493, \ 0.1752, \ 0.1027, 
              0.1124}}, 
blocks\  = \ 
          StringPartition[s, \ m]; 
counts\  = \ 
          Map[\ Max[1, \ 
                Min[4, Max[StringLength[StringSplit[#, \ "0"]]]]] &, \ 
            blocks]; 
chisqr\  = \ 
          Apply[Plus, 
            MapThread[\ \((\((#1\  - \ 
                          n*\ #2)\)\ ^2)\)\ /\ \((n\ *\ #2)\) &, \ {counts, \ 
                pik}]]; 
pval\  = \ 
          GammaRegularized[k/2, \ 
            chisqr\ /\ 2]; 
pval\  ≥ \ 
          0.01
]


    LongestRunsOnes10000[
        s_String]\  := 
Module[{blocks, m = 10^4, \ 
          k = 6, \ n = \ 75, \ 
          pik\  = \ {0.0882, \ 0.2092, \ 0.2483, \ 0.1933, \ 0.1208, 
              0.0675, \ 0.0727}}, 
blocks\  = \ 
          StringPartition[s, \ m]; 
counts\  = \ 
          Map[\ Max[1, \ 
                Min[4, Max[StringLength[StringSplit[#, \ "0"]]]]] &, \ 
            blocks]; 
chisqr\  = \ 
          Apply[Plus, 
            MapThread[\ \((\((#1\  - \ 
                          n*\ #2)\)\ ^2)\)\ /\ \((n\ *\ #2)\) &, \ {counts, \ 
                pik}]]; 
pval\  = \ 
          GammaRegularized[k/2, \ 
            chisqr\ /\ 2]; 
pval\  ≥ \ 
          0.01
]\)\)\)

2.5 Binary Matrix Rank Test

(*2.5*) 
BinaryMatrixRankTest[
          s_String,  m_:  32,  q_:  32]  :=  
        Module[{}, 
r  = 
            Map[MatrixRank[#,  Modulus -> 2] &, 
              Partition[Partition[(ToCharacterCode[s]  -  48),  m],  
                q]]; 
n  =  
            Length[r]; 
fm  =  
            Count[r,  m]; 
fm1  =  
            Count[r,  
              m - 1]; 
chisqr  =  ((fm  -  
                        0.2888 n)^2)/(0.2888  
                    n)  +  ((fm1  -  
                        0.5776 n)^2)/(0.5776 n)  +  ((n  -  
                        fm  - fm1  -  
                        0.1336  
                          n)^2)/(0.1336 n); 
pval  =  Power[E,  -chisqr / 
                2]; 
pval  ≥  
            0.01
];

2.6 Spectral Test

The spectral test can be visualized after running the test using

ListPlot[dstreal, Axes->False, Frame->True, PlotRange->{Automatic, {t, Max[dstreal]}}]

(*2.6*) 
SpectralTest[str_String]  :=  
        Module[{n  =  StringLength[str],  str2,  
            dst}, 
str2  =  (ToCharacterCode[str]  -  
                    48)*2 - 1; 
dst  =  
            Fourier[str2,   FourierParameters -> {1, -1}]; 
dstreal  = Take[Abs[dst],  n/2]; 
t = 
            Sqrt[3 *n]; 
n0  =  
            0.95 *n / 2; 
n1  =  
            Length[Select[
                dstreal,  # < t &]]; 
d =  (n1  -  
                  n0) / 
              Sqrt[n * 0.95 * 0.05 / 2]; 
pval  =  
            Erfc[Abs[d] / 
                Sqrt[2]]; 
pval  ≥  
            0.01
];

2.7 Non-overalpping Template Matching Test

(*2.7*) 
 (* 
s  =  
        input, 
b  =  template string
M  =  
            blocks length
*) 
NonOverlappingTemplateMatching[s_String,  b_String,  
          M_:  131072]  := 
Module[{m  =  
              StringLength[
                b]}, 
 (* if b  <  2  ||  
                  b  >  
                    21) error *) 
blocks  =  
            StringPartition[s, M]; 
counts  =  
            Map[StringCount[#,  b] &,  
              blocks]; 
avg  =  
            N[(M  -  m + 1)/(2^m)]; 
var  =  
            N[M*(2^-m  -  (2*m  - 
                          1)*2^(-2*
                            m))]; 
chisqr  =  
            Sum[(counts[[i]]  -   avg)^2,  {i, 1,  
                  Length[counts]}] / var; 
pval  =  
            GammaRegularized[Length[blocks]/2,  
              chisqr / 2]; 
pval  ≥  
            0.01
];

2.8 Overlapping Template Matching Test

\(\(\( (*2.8*) \)\(
\)\(
\)\( (*\ 
      crazy\ function\ to\ compute\ values\ *) \)\(
\)\(Pr[
        0, \ x_]\  := \ Power[E, \ \(-x\)]

    \(Pr[u_Integer, \ x_]\  := 
        x*\((E^\(-x\))\)*\((2^\(-u\))\)*\ 
          Sum[\ Binomial[u - 1, \ i\  - \ 1]*\((x^\((i - 1)\))\)/
                Factorial[i], \ {i, \ 1, 
              u}];\)
 (*\(\(\(\ \)\(OR

          
          Pr[u_Integer, \ 
              x_]\  := \ \ x*\((E^\(-\((2  x)\)\))\)*\((2^\(-u\))\)*
              Hypergeometric1F1[u + 1, 2, \ 
                x]\)\(]\)\);\)
*) 


    \(OverlapingTemplateMatching[s_String, \ b_String, \ M_:  1032, \ 
          K_:  5]\  := \ 
Module[{}, \ 
n\  = \ StringLength[s]; 
bigN\  = \ 
            Floor[n/\ M]; 
m\  = \ 
            StringLength[b]; 
lambda\  = \ 
            N[\((M - \ \ m\  + 
                    1)\)*\((2^\(-m\))\)]; 
eta\  = \ 
            lambda\ /\ 2; 
pi\  = \ 
            Table[Pr[i, \ eta], {i, \ 0, \ 
                K - 1}]; 
AppendTo[pi, \ 
            1 - \ Apply[Plus, \ pi]]; 
v\  = \ 
            Table[0, \ {K + 
                  1}]; 
 (*\(\(\ \)\(if\ b\  < \ 2\  || \ 
                  b\  > \ 
                    21\)\()\)\)\ error\ *) 
blocks\  = \ 
            StringPartition[s, \ M]; 
blockLen\  = \ 
            StringLength[blocks[\([1]\)]]; 
 (*\ 
            count\ occurances, \ but\ if\ more\ than\ K, \ 
            just\ make\ it\ K\ *) 
counts\  = \ 
            Map[Min[StringCount[#, \ b, \ Overlaps -> True], \ K] &, \ 
              blocks]; 
Scan[\((\(v[\([# + 
                        1]\)]++\))\) &, \ 
            counts]; 
chisqr\  = \ 
            Sum[\ \((\((v[\([i]\)]\  - \ bigN*\ pi[\([i]\)]\ )\)^2)\)/\((bigN*
                    pi[\([i]\)])\), \ \ {i, 1, 
                K + 1}]; 
pval\  = \ 
            GammaRegularized[K/2, \ 
              chisqr\ /\ 2]; 
pval\  ≥ \ 
            0.01
];\)\)\)\)

2.9 Marurer's Universal Statistic Test

(*2.9*) 
MaurersUniversalStatisticTest[
        s_String,  L_Integer,  
        q_Integer]  := 
Module[
 (* 
          Table comes from  section  5.4  .5 of Handbook of Applied Cryptography *) {ru  =  {{0.7326495,  0.690},  {1.5374383,  
                1.338},  {2.4016068,  1.901},  {3.3112247,  
                2.358},  {4.2534266,  2.705},  {5.2177052,  
                2.954},  {6.1962507,  3.125},  {7.1836656,  
                3.238},  {8.1764248,  3.311},  {9.1723243,  
                3.356},  {10.170032,  3.384},  {11.168765,  
                3.401},  {12.168070,  3.410},  {13.167693,  
                3.416},  {14.167488,  3.419},  {15.167379,  
                3.421}}}, 
 (* 
          TURN INTO LIST OF NUMBERS 1  -  
            2^L *) 
blocks  =  
          Map[1  +  FromDigits[#,  2] &,  
            Partition[ToCharacterCode[s]  -  48,  
              L]]; 
k  =  
          Length[blocks]  -  q; 
states  =  
          Table[0,  {i, 2^L}]; 
 (* 
          do initialization *) 
Do[ 
          states[[ blocks[[i]] ]]  =  i ,  {i, 1,  
            q}]; 
sum =  
          0.0; 
Do[  sum  += 
            Log[2,  i  -  
                states[[ 
                    blocks[[i]] ]]]; 
states[[ 
                blocks[[i]] ]]  =  i,  {i,  q + 1,  
            Length[blocks]}
]; 
fn  =  
          sum / k; 
pval  =  
          Erfc[Abs[fn  -  ru[[L,  1]] ] / 
              Sqrt[2 * 
                  ru[[L, 
                      2]] ] ]; 
pval  ≥  0.01
]

2.10 Lempel-Ziv Compression Test

\(\(\( (*2.10*) \)\(
\)\(
\)\( (*\ 
      This\ is\ the\ slow\ but\ simple\ way . \ 
          Good\ for\ illustration\ *) \)\(
\)\(LempelZivCompressionTestSlow[
          s_String]\  := 
Module[{i\  = \ 1, \ j\  = \ 0, \ 
            n\  = \ StringLength[s], \ tmp, \ μ\  = \ 69586.25, \ 
            sigma\  = \ 
              70.448718}, 
words\  = \ {}; 
While[\ 
            i\  + j\ \  ≤ \ n, 
tmp\  = \ 
              StringTake[s, \ {i, \ i + j}]; 
If[
              MemberQ[words, \ 
                tmp], \ 
\(j++\), 
\((AppendTo[words, \ tmp]; \ i\  += \ j + 1; \ 
                j\  = \ 
                  0)\)];
]; 
wobs\  = \ 
            Length[words]; 
pval\  = \ 
            Erfc[\((μ\  - \ wobs\ )\)/Sqrt[2*sigma]]\ /
              2; 
pval\  ≥ \ 
            0.01\ 
];\)\)\)

Instead of linear searching, this uses the "hash-table" trick,
http://library.wolfram.com/conferences/devconf99/lichtblau/Links/index_lnk_5.html

LempelZivCompressionTest[
          s_String]  := 
Module[{minj  =  0, 
            i  =  1,  j  =  0,  n  =  StringLength[s],  
            tmp,  μ  =  69586.25,  sigma  =  70.448718,  words,  
            wordlen}, 
s2  =  
            BinaryStringToDigits[s]; 
wobs  =  
            0; 
While[ 
            i  + j   ≤  
              n, 
tmp  =   s2[[ 
                  Range[ i,  i + j] ]]; 
If[
              Head[words[tmp]]  =!=  words,  j++;  
              Continue[]]; 
words[tmp]  =  
              tmp; 
 ++wobs; 
i  +=  j + 1; 
 (* optimization,  
              if we got all the words of a certain length, 
then skip over them .   This shaves 30  % off the runtime *) 
If[Head[wordlen[j]]  ===  wordlen, 
              wordlen[j]  = 
                Power[2,  j + 1] - 1,  --wordlen[
                  j]]; 
If[
              wordlen[minj]  ==  
                0,  ++minj]; 
 
            j  =  minj
]; 
pval  =  
            Erfc[(μ  -  wobs )/Sqrt[2*sigma]] /
              2; 
pval  ≥  
            0.01 
];

2.11 Linear Complexity Test

\(\(\( (*2.11*) \)\(
\)\(
\)\(
\)\( (*\ Just\ computes\ the\ linear\ complexity . 
\ \ Compiling\ this\ gives\ 10  
        x\ boost\ *) \)\(
\)\(\(LinearComplexity\  = \ 
        Compile[{{u, \ _Integer, \ 1}}, \ 
          Module[{len\  = \ Length[u], \ b, c, d, p, tmp, \ l = 0, \ 
              m = 0, \ }, 
c = 
              ReplacePart[Table[0, \ {len}], \ 1, 1]; 
b = 
              ReplacePart[Table[0, \ {len}], \ 1, 
                1]; 
Do[
d = \ 
                Mod[u[\([n]\)]\  + \ 
                    Dot[\ Take[c, \ {2, \ l + 1}], 
                      Reverse[Take[u, {n - l\ , \ n - 1\ }]]], \ 
                  2]; 
If[
                d == 1, 
tmp\  = \ 
                  c; 
p\  = \ 
                  Table[0, \ {len}]; 
Do[\ \ If[
                    b[\([i]\)]\  == 1, 
                    p[\([i + n - m]\)]\  = \ 1], \ {i, 1, \ 
                    l + 1}]; 
c\ \  = \ 
                  Mod[c + p, \ 2]; 
If[
                  2*l\  ≤ \ n, 
l\  = \ 
                    n\  - l; 
m\  = \ 
                    n; 
b\  = \ 
                    tmp;
];
], \ {n, \ 
                1, \ len}]; 
l
]
];\)

 (*\ 
      m\  = \ length\ of\ bits\ in\ block\ *) 

    \(LinearComplexityTest[s_String, \ 
          m_Integer]\  := 
Module[{k = 6, \ 
            pi = {0.01047, \ 0.03125, \ 0.125, 0.5, 0.25, 0.0625, \ 
                0.02078}}, 
 (*\ 
            theoretical\ average\ of\ what\ the\ Linear\ Complexity\ *) 
avg\  = \ \((m/
                  2)\)\  + \ \((9\  + \ Power[\(-1\), \ m\  + 1])\)/
                36\  - \ \((m/3\  + \ 2/9)\)/
                Power[2, m]; 
 (*\ 
            blocks\ *) 
blocks\  = \ 
            StringPartition[s, m]; 
 (*\ 
            number\ of\ blocks\ *) 
bigN\  = \ 
            Length[blocks]; 

 (*\ 
            get\ lc\ for\ each\ block\ *) 
lc\  = \ 
            Map[LinearComplexity[BinaryStringToDigits[#]] &, \ 
              blocks]; 
chisqr\  = \ 
            Sum[Power[v[\([i]\)]\  - \ bigN\ *pi[\([i]\)], 
                  2]\ /\ \((bigN\ *\ pi[\([i]\)])\), \ {i, 1, \ 
                Length[v]}]; 
t\  = \ 
            Map[\((Power[\(-1\), \ m]*\((#\  - \ avg)\)\  + \ 2/9)\) &, \ 
              lc]; 
 (*\ 
            spec\ says\ Ti\  ≤ \ \(-2.5\), \ etc, \ 
            while\ RangeCount\ is\ ti\  < \ \(-2.5\), \ 
            should\ make\ no\ practical\ difference\ *) 
v\  = \ Reverse[\ 
              RangeCounts[\(-1\)*t, \ {\ \(-2.5\), \ \(-1.5\), \ \(-0.5\), \ 
                  0.5, \ 1.5, \ 2.5}]]; 
chisqr\  = \ 
            Sum[\ Power[\ v[\([i]\)]\  - \ bigN*pi[\([i]\)], 
                  2]/\ \((bigN\ *\ pi[\([i]\)])\), \ {i, 1, \ 
                Length[v]}]; 
pval\  = \ 
            GammaRegularized[k/2, \ 
              chisqr\ /\ 2]; 
pval\  ≥ \ 
            0.01
];\)\)\)\)

2.12 Serial Test

\(\(\( (*2.12*) \)\(
\)\( (*\ 
      TODO : \ change\ psmi\ to\ something\ more\ natural\ looking . \ It' 
          s\ the\ sum\ of\ squares\ of\ the\ first\ part\ of\ the\ lst\ *) \)\(
\)\( (*\ 
      TODO : \ make\ stringpartition\ do\ the\ wrap\ around\ *) \)\(
\)\(SerialTest[s_String, \ 
        m_Integer]\  := 
Module[{}, 
n\  = \ StringLength[s]; 
f1\  = \ 
          Frequencies[
            StringPartition[StringJoin[s, StringTake[s, m - 1]], m, 
              1]]; 
f2\  = \ 
          Frequencies[
            StringPartition[StringJoin[s, \ StringTake[s, m - 2]], m - 1, 
              1]]; 
f3\  = \ 
          Frequencies[
            StringPartition[StringJoin[s, \ StringTake[s, m - 3]], m - 2, 
              1]]; 
psim1\  = \ 
          0; 
psim2\  = \ 0; 
psim3\  = \ 
          0; 

If[
          m\  ≥ \ 0, 
          psim1\  = \ 
            Power[2, 
                  m]\ *\ \ Fold[Plus[#1, \ Power[Part[#2, 1], \ 2]] &, \ 0, \ 
                    f1]\ /n\  - n]; 
If[
          m\  ≥ \ 1, 
psim2\  = \ 
            Power[2, 
                  m - 1]\ *\ \ Fold[Plus[#1, \ Power[Part[#2, 1], \ 2]] &, \ 
                    0, \ f2]\ /n\  - n]; 
If[
          m\  ≥ \ 2, 
psim3\  = \ 
            Power[2, 
                  m - 2]\ *\ \ Fold[Plus[#1, \ Power[Part[#2, 1], \ 2]] &, \ 
                    0, \ f3]\ /n\  - n]; 
d1\  = \ 
          psim1\  - \ psim2; 
d2\  = \ 
          psim1\  - \ 2*psim2\  + \ 
            psim3; 
pval\  = \ {GammaRegularized[
              2^\((m - 2)\), d1/2], \ \ GammaRegularized[2^\((m - 3)\), 
              d2/2]}; 
Map[\((#\  ≥ \ 
                0.01)\) &, \ pval]
]\)\)\)

2.13 Approximate Entropy Test

ApproximateEntropyTest[
        s_String,  
        m_Integer]  := 
Module[{n  =  
            StringLength[s]}, 
f1  =  
          Frequencies[
            StringPartition[StringJoin[s, StringTake[s, m - 1]], m, 
              1]]; 
f2  =  
          Frequencies[
            StringPartition[StringJoin[s,  StringTake[s, m]], m + 1, 
              1]]; 
c1  =  
          Map[#[[1]] &,  f1] / n; 
c2  =  
          Map[#[[1]] &,  f2] / n; 
phi1  =  
          Apply[Plus,  
            Map[# Log[#] &,  c1]]; 
phi2  =  
          Apply[Plus,  
            Map[# Log[#] &,  c2]]; 
apen  =  
          phi1  -  phi2; 
chisqr  =  
          2  n (Log[2]  - apen); 
pval  = 
          GammaRegularized[2^(m - 1), 
            chisqr/2]; 
pval  ≥  
          0.01
]

2.14 Cumulative Sums Test

CumulativeSumsTest[
        s_String]  :=  
      Module[ {n  =  StringLength[s], cs,  
          ndist  =  NormalDistribution[0, 1]}, 
cs  =  
          CumulativeSums[ (ToCharacterCode[s]  -  48)*2 - 
              1]; 
z  = 
          Max[Abs[cs]]; 
pval =  
          1 - Apply[Plus,  
              Table[CDF[ndist,  (4  k + 1) z/Sqrt[n]]  -  
                  CDF[ndist,  (4  k - 1) z/Sqrt[n]],  {k,  
                  Floor[(-Floor[n/z]  +  1)/4],  
                  Floor[(Floor[n/z]  - 1)/4]}]]  +  
            Apply[Plus, 
              Table[CDF[ndist, (4  k + 3) z/Sqrt[n]]  -  
                  CDF[ndist,  (4  k + 1) z/Sqrt[n]],  {k,  
                  Floor[(-Floor[n/z]  - 3)/4],  
                  Floor[(Floor[n/z]  - 1)/
                      4]}]]; 
pval  ≥  
          0.01
]

2.15 Random Excursions Test

\(\(\( (*2.15*) \)\(
\)\( (*\ given\ a\ list\ of\ cycles, \ 
      count\ the\ number\ of\ times\ x\ appears\ in\ each\ *) \)\(
\)\(\(freqCycles[cycles_List, \ x_Integer]\  := \ 
        Module[{tmp}, \ 
          tmp\  = \ 
            Map[Min[{5, \ Count[#, \ x]}] &, \ cycles]; 
\ 
          Table[Count[tmp, \ i], \ {i, 0, 5}]
];\)

    \(Remove[pik];\)

    \(pik[0, \ x_]\  := \ 1\  - \ 1/\((2*Abs[x])\);\)

    \(pik[k_, \ 
          x_]\  := \ \((1/\((4*x*x)\))\) \((1 - 1/\((2*Abs[x])\))\)^\((k - 
                1)\);\)

    \(pik[5, \ 
          x_]\  := \ \((1/\((2*Abs[x])\))\) \((1 - 1/\((2*Abs[x])\))\)^4;\)

    \(xvals\  = \ {\(-4\), \(-3\), \(-2\), \(-1\), 1, 2, 3, 4};\)

    \(pikTable\  = \ Map[Table[\ pik[i, \ #], \ {i, 0, 5}] &\ , \ xvals];\)

    \(achisqr[x_, \ y_, \ j_]\  := \ 
        Apply[Plus, 
          MapThread[\ \((\((#1\  - \ j*#2)\)^2)\)/\((j\ *\ #2)\) &, \ {x, \ 
              y}]];\)

    \(RandomExcursionsTest[
          str_String]\  := 
Module[{}, \ 
 (*\ make\ cumulative\ sums, \ 
            and\ add\ zeros\ at\ start\ and\ end\ *) 
cumsum\  = \ 
            Flatten[{0, 
                CumulativeSums[\((ToCharacterCode[str]\  - \ 48)\)*2 - 
                    1\ ], \ 0}]; 
 (*\ 
            Find\ where\ the\ zeros\ are, \ 
            and\ group . 
\ The\  + 
              1, \ \(-1\)\ are\ the\ positions\ of\ start\ and\ end\ of\ a\ cycle*) 
cyclePositions\  = \ 
            Partition[Flatten[Position[cumsum, \ 0]], \ 2, \ 
              1]; 
 (*\ 
            make\ a\ list\ of\ cycles\ *) 
cycles\  = \ 
            Map[Take[cumsum, \ #] &, \ 
              cyclePositions]; 
j\  = \ 
            Length[cycles]; 
stateCyclesTable\  = \ 
            Map[freqCycles[cycles, #] &, \ xvals]; 
 (*\ 
            with\ the\ freq\ cycles, \ 
            and\ the\ pi\ table\ compute\ chisqr\ for\ a\ state*) 
chisqr\  = 
            MapThread[
              achisqr[#1, \ #2, \ j] &, \ {stateCyclesTable, \ 
                pikTable}]; 
pval\  = \ 
            Map[GammaRegularized[\ 5/2, \ #/2] &, \ 
              chisqr]; 
Map[#\  ≥ \ 0.01 &, \ 
            pval]
];\)\)\)\)

2.16 Random Excusions Variant Test

(*2.16*) 
 (* 
      helper function *) 
getFreq[f_List,  
          x_]  :=  Module[{result}, 
result  =  
            Select[f,  #[[2]]  ==  x  &]; 
If[
            Length[result] == 0,  0,  
            result[[1, 
                1]]]
];


    RandomExcursionsVariantTest[
        s_String]  := 
Module[{},  
        cumsum  =  
          CumulativeSums[(ToCharacterCode[s]  -  48)*2 - 
              1 ]; 
freq  =  
          Select[Frequencies[cumsum],  
            Abs[#[[2]] ]  ≤  
                9  &]; 
j  =  
          getFreq[freq,  0] + 1; 
pval  = 
          Map[ {#, 
                Erfc[Abs[getFreq[freq, #]  -  j]/
                    Sqrt[2 *
                        j*(4*Abs[#]  -  
                            2)]]} &,  {-9, -8, -7, -6, -5, -4, -3, -2, -1, 1, 2, 3, 4, 5, 6, 7, 8, 
              9}]; 
Map[{#[[1]], #[[2]] ≥  0.01} &,  pval]
]

Change Log

26-Jun-05

Initial release 0.90

06-Jul-05 - 0.9.1

> 20x improvement in LinearComplexity!  The test should now run in under 3 minutes. The symbolic computations were replaced with custom (more c-like) arrays.  This provided 2x.  Then compiling the function added 10x more!

Cleanup of some variable names

Serial Test cleaned up to return both P-values as a list

07-Oct-06 0.9.2

Cleaned up Spectral Test and a few others.