In [1]:
import numpy as np
from scipy.stats import chisqprob as cdchisq
from fgn import fbm
from xtree import *
import matplotlib.pyplot as plt
%matplotlib inline
In [2]:
generator = fbm( N = 2**20+1, H = .65 )
generator.set_rnd( np.random.RandomState( 123 ) )
T, X = generator( )
In [3]:
delta = np.std( np.diff( X ) )
In [4]:
hp = list( ) ; ht = list( ); hx = list( ) ; ex = list( ) ; wt = list( )
delta = np.std( np.diff( X ) ) if delta is None else delta
Z = ( X - X[ 0 ] ) / delta
lht0, lhp0, lhx0, lex0, lwt0 = xtree_integer_crossings( T, Z )
height = 1 ; delta *= 2
Z = ( X - X[ 0 ] ) / delta
lht1, lhp1, lhx1, lex1, lwt1 = xtree_integer_crossings( T, Z )
In [51]:
hits = np.searchsorted( lht0, lht1 )
last_hit = hits // 2 - 1
directions = np.sign( np.diff( lhp0 ) )[ 1::2 ]
up_down_mask = directions < 0
up_down_mask[ last_hit[ 1: ] ] = False
up_down_total = np.cumsum( up_down_mask )[ last_hit ] ; up_down_total[ 0 ] = 0
excursions = np.empty( ( len( hits ) - 1, 3 ), np.int )
excursions[:,0] = np.diff( up_down_total )
excursions[:,1] = np.diff( hits ) // 2 - excursions[:,0] - 1
excursions[:,2] = directions[ last_hit[ 1: ] ]
In [60]:
b = list( )
ud = 0 ; du = 0
directions = np.sign( np.diff( lhp0 ) )
ds = "".join( "-+"[int(d+1)//2] for d in directions )
for i in range( 0, len( ds ), 2 ) :
pp = ds[i:i+2]
if pp == '+-' :
ud += 1
elif pp == '-+' :
du += 1
else :
b.append( ( ud, du, pp ) )
ud = du = 0
In [66]:
a = list( )
ud = 0 ; du = 0
directions = np.sign( np.diff( lhp0 ) )
for f,s in zip( directions[0::2], directions[1::2] ) :
if f > 0 and s < 0 :
ud += 1
elif f < 0 and s > 0 :
du += 1
else :
a.append( ( ud, du, s ) )
ud = du = 0
aa = np.array( a )
In [68]:
print np.allclose( aa, excursions )
In [63]:
print ds[:26]
print b[:10], a[:10]
print excursions[:10]
print excursions.shape
In [31]:
uu = np.sum( ( directions[1::2] > 0 ) & ( directions[0::2] > 0) )
du = np.sum( ( directions[1::2] > 0 ) & ( directions[0::2] < 0) )
ud = np.sum( ( directions[1::2] < 0 ) & ( directions[0::2] > 0) )
dd = np.sum( ( directions[1::2] < 0 ) & ( directions[0::2] < 0) )
In [33]:
print excursions.sum( axis = 0)
In [36]:
print uu, ud, du, dd
print uu+ud, u_
print uu+du, _u
print du+dd, d_
print ud+dd, _d
In [30]:
print _d, _d-dd, d_, dd
print _u, _u-uu, u_, uu
In [ ]:
max_levels, max_crossings = 10, 10
delta = np.std( np.diff( X ) )
Tnk, hp, Znk, ex, Wt = xtree_build( T, X, delta = delta )
# Tnk, hp, Znk, ex, Wt = xtree_build( T, X, delta = delta, max_height = max_levels + 1 )
In [ ]:
In [ ]:
counts = np.zeros( ( max_levels + 1, max_crossings // 2 ), dtype = np.float )
for level, subcrossings in enumerate( Znk[ 1: ], 0 ) :
sx_count, sx_freq = np.unique( subcrossings, return_counts = True )
sx_count = np.minimum( sx_count, max_crossings )
level = np.minimum( level, max_levels )
counts[ level, sx_count // 2 - 1 ] = sx_freq
In [ ]:
Nn = np.array( [len( Tk ) for Tk in Tnk ], dtype = np.int ).reshape( ( len( Tnk ), -1 ) )
In [ ]:
levels = range( 1, np.minimum( 4, len( Nn ) ) )
p_count = np.sum( counts[ levels,: ], axis = 0 )
p_distr = p_count / p_count.sum( )
distr = counts / counts.sum( axis = 1 ).reshape( counts.shape[ 0 ], -1 )
stat = np.sum( Nn[ levels ] * ( ( distr[ levels,: ] - p_distr )**2 / p_distr ) )
dof = len( levels ) * ( distr.shape[ 1 ] - 1 )
cdchisq( stat, dof )
In [ ]:
[ cdchisq( f/2.0, f ) for f in range( 1,10 )]
In [ ]: